Take-home Exercise 1 - Part 1

Harnessing Geospatial Analytics to Uncover Armed Conflict Patterns in Myanmar

Author

Foo Jia Yi Samantha

Published

September 22, 2024

Modified

September 22, 2024

1. Overview

1.1 Background of Myanmar’s Long-Standing Conflicts

The conflict in Myanmar is not just a result of the coup but is deeply rooted in the country’s decades-old complex ethnic and political landscape, characterised by tensions between the central government and various ethnic minority groups, each with its own armed forces. The post-coup violence has exacerbated these long-standing conflicts, leading to a severe humanitarian crisis, with thousands killed, hundreds of thousands displaced, and widespread human rights abuses reported.

1.2 Objectives of Take-home Exercise 1

As such, Geospatial analytics has become a valuable tool for evaluating and comprehending the intricacies of increasing conflicts. This exercise aims to reveal the spatial and spatio-temporal distribution of armed conflict in Myanmar by leveraging spatial point pattern analysis. Additionally, it aims to gain clearer insights into the geographical and logistical patterns of violence throughout the nation.

By the end of this take-home exercise, I aim to complete these steps in my spatial point pattern analysis in uncovering the distribution of armed conflict in Myanmar.

  • Using appropriate function of sf and tidyverse packages, import and transform the downloaded armed conflict data and administrative boundary data into sf tibble data.frames.
  • Using the geospatial data sets prepared, derive quarterly KDE layers.
  • Using the geospatial data sets prepared, perform 2nd-Order Spatial Point Patterns Analysis.
  • Using the geospatial data sets prepared, derive quarterly spatio-temporal KDE layers.
  • Using the geospatial data sets prepared, perform 2nd-Order Spatio-temporal Point Patterns Analysis.
  • Using appropriate tmap functions, display the KDE and Spatio-temporal KDE layers on openstreetmap of Myanmar.
  • Describe the spatial patterns revealed by the KDE and Spatio-temporal KDE maps.

1.3 About the Datasets

1) Armed Conflict Data (From ACLED)

This Armed Conflict Location & Event Data (ACLED) is an independent, impartial, international non-profit organisation which owns an extensive database of violent conflict and protest in countries and territories around the world.

For the purpose of this exercise, I have downloaded ACLED’s data on Myanmar which includes a series of conflict events, particularly between 1 January 2021 to 30 June 2024.

🔗 Source: ACLED

📁 Format: comma separated values (CSV)

As the dataset is rather extensive, I will be performing my analysis on armed conflict events in a quarterly basis to streamline my tasks. The data included in this dataset are as follows:

Event Type

ACLED categorises events into various types. I will mainly be focusing on these four event types: Battles, Explosion/Remote violence, Strategic developments, and Violence against civilians.

  • event_id_cnty: unique ID for each conflict
  • event_type: category of event e.g. Battle, Violence Against Civilians, Protests, Explosions/Remote Violence, Strategic Developments
  • sub_event_type: a more detailed classification within event type
  • disorder_type: classifies the event based on the nature of the disorder e.g. political violence, demonstrations, strategic developments[A1]
  • civilian_targeting: yes/no value, whether event involves specifically targeting civilians
  • Note: when “strategic developments” are used in Event Type, it is also used in the disorder type (vice-versa)

Location and Geospatial Data

The database provides detailed geographic information, pinpointing the exact or approximate locations of conflict events across Myanmar. This includes cities, towns, and rural areas.

  • iso: the country code for Myanmar which uses 104 in this case
  • region: region of conflict within Myanmar
  • country: indicates Myanmar
  • admin1, admin2, admin3: 1st, 2nd and 3rd level administration division within Myanmar e.g. states, division, sub-division
  • location: specific geographic location or name of the place where the conflict event occurred
  • latitude: latitude of the conflict event
  • longitude: longitude of the conflict event
  • geo_precision: indicates the level of precision for the geographic coordinates provided

Date and Time

ACLED records the specific dates and, where possible, times of conflict events.

  • event_date: date of conflict
  • year: year of conflict
  • time_precision: accuracy of the date and time information provided

Actors

  • Indicate the actors involved in the conflict, such as the Tatmadaw (Myanmar’s military), ethnic armed organizations, local militias, civilian protestors, and other groups.
  • actor1: primary actor involved in the conflict event. E.g. a government force, rebel group, militia, or any organised entity
  • assoc_actor_1: a secondary group that is aligned with or supports the primary actor (Actor1) in the event
  • inter1: an interaction code that categorises actor1, could be a government force, rebel group, military force, rioter, civilian, or other entities
  • interaction: combined description of actor1 and actor2 (no particular order of aggression)

Fatalities

  • fatalities: tracks the number of reported fatalities associated with each conflict event

Others

  • source: source of information for the conflict event
  • source_scale: scale of the source e.g. local, national, international
  • notes : additional comments
  • tags: keywords associated with the conflict event
  • timestamp: date and time when conflict event was entered/updated in the database

2) Geospatial Data (From Myanmar Information Management Unit)

I will also be using a geospatial dataset from the Myanmar Information Management Unit (MIMU) in shapefile (.shp) format, specifically of the Myanmar state at the 2nd administrative level with district boundaries.

🔗 Source: MIMU

📁 Format: shapefile (.shp)

My reason for choosing the district boundary dataset is that we do not want to select a boundary dataset that is too generalised when analysing conflict events since it might not provide sufficient insights to trends where conflict events happen. Neither do we want to analyse a geography that is too detailed (e.g. Admin 3 - townships) since it can be computationally inefficient as seen in the types of boundary data below.

I have donwloaded the two data sets and organised them into my folder as follows.

2. Let’s Set Up!

2.1 Importing Libraries into R

To carry out this exercise, I will be using the following R packages:

  • sf: a relatively new R package specially designed to import, manage and process vector-based geospatial data in R.
  • spatstat: has a wide range of useful functions for point pattern analysis. In this take-home exercise, it will be used to perform 1st- and 2nd-order spatial point patterns analysis and derive kernel density estimation (KDE) layer.
  • raster: reads, writes, manipulates, analyses and model of gridded spatial data (i.e. raster). In this take-home exercise, it will be used to convert image output generate by spatstat into raster format.
  • tmap / sparr: provides functions for plotting cartographic quality static point patterns maps or interactive maps by using leaflet API.
  • tidyverse / dplyr: for transforming and organising data for analysis
  • magick: used for plotting animated map plots into GIFs for the spatio-temporal point pattern analysis

Now, let’s install and load these packages in RStudio.

pacman::p_load(sf, raster, spatstat, sparr, tmap, tidyverse, magick, dplyr)

2.2 Importing Data Sets into R

1) Armed Conflicts Data

Next, I will import the downloaded armed conflict data. For aspatial datasets like this, we will import into Rstudio using read_csv() function of the readr package.

# Import armed conflict data
conflict_data <- read_csv("data/aspatial/2021-01-01-2024-06-30-Myanmar.csv")
Rows: 87746 Columns: 28
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr (18): event_id_cnty, event_date, disorder_type, event_type, sub_event_ty...
dbl (10): year, time_precision, inter1, interaction, iso, latitude, longitud...

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Observations

The 2021-01-01-2024-06-30-Myanmar.csv dataset contains 87746 rows and 28 columns which indicates the presence of 87746 unique armed conflict events in Myanmar.

After importing the dataset, we can inspect the dataset using the glimpse() function.

# Inspect the conflict data
glimpse(conflict_data)
Rows: 87,746
Columns: 28
$ event_id_cnty      <chr> "MMR64313", "MMR64313", "MMR64320", "MMR64320", "MM…
$ event_date         <chr> "30 June 2024", "30 June 2024", "30 June 2024", "30…
$ year               <dbl> 2024, 2024, 2024, 2024, 2024, 2024, 2024, 2024, 202…
$ time_precision     <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, …
$ disorder_type      <chr> "Political violence", "Political violence", "Politi…
$ event_type         <chr> "Battles", "Battles", "Battles", "Battles", "Battle…
$ sub_event_type     <chr> "Armed clash", "Armed clash", "Armed clash", "Armed…
$ actor1             <chr> "People's Defense Force - Mandalay", "Military Forc…
$ assoc_actor_1      <chr> "MDA - AGF: Madaya - The Authentic Genes Force; SST…
$ inter1             <dbl> 3, 1, 3, 1, 3, 1, 1, 3, 1, 1, 1, 2, 2, 1, 1, 2, 1, …
$ interaction        <dbl> 13, 13, 13, 13, 13, 13, 10, 13, 13, 10, 12, 12, 12,…
$ civilian_targeting <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,…
$ iso                <dbl> 104, 104, 104, 104, 104, 104, 104, 104, 104, 104, 1…
$ region             <chr> "Southeast Asia", "Southeast Asia", "Southeast Asia…
$ country            <chr> "Myanmar", "Myanmar", "Myanmar", "Myanmar", "Myanma…
$ admin1             <chr> "Mandalay", "Mandalay", "Mandalay", "Mandalay", "Ma…
$ admin2             <chr> "Mandalay", "Mandalay", "Pyinoolwin", "Pyinoolwin",…
$ admin3             <chr> "Patheingyi", "Patheingyi", "Singu", "Singu", "Thab…
$ location           <chr> "Aung Tha Pyay", "Aung Tha Pyay", "Pin Lel Gyi", "P…
$ latitude           <dbl> 22.1504, 22.1504, 22.5752, 22.5752, 22.8800, 22.880…
$ longitude          <dbl> 96.2364, 96.2364, 96.0661, 96.0661, 95.9700, 95.970…
$ geo_precision      <dbl> 2, 2, 1, 1, 2, 2, 2, 2, 2, 2, 1, 1, 1, 1, 2, 2, 1, …
$ source             <chr> "Democratic Voice of Burma; Irrawaddy", "Democratic…
$ source_scale       <chr> "National", "National", "National", "National", "Na…
$ notes              <chr> "On 30 June 2024, near Aung Tha Pyay village (Pathe…
$ fatalities         <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 2, 0, 0, 0, …
$ tags               <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,…
$ timestamp          <dbl> 1720552468, 1720552468, 1720552468, 1720552468, 172…
Observations

The event_date field shows that it uses a character datatype instead of date - we will fix this later. Also, we can observe that thelongitude and langitude fields appear to be adopting the WGS84 geographic coordinate system since they are in the -180/180 and -90/90 range respectively.

2) Myanmar Boundary Data

Observations

When working with Myanmar’s boundary, we need to assign the appropriate coordinate reference system. However, since Myanmar is split into two UTM - West Myanmar (crs: 32646) and East Myanmar (crs: 32647).

Hence, I will also import the administrative boundary data into a simple features tibble data.frame using st_read() of the sf package and check the number of rows returned for both CRS 32646 and 32647. This function reads the shapefile data and returns an sf object that can be used for further analysis.

Find out conflicts count by CRS
conflict_crs <- st_as_sf(conflict_data, coords = c("longitude", "latitude"), crs = 4326) 

# Count number of conflicts for CRS 32646
conflict_data_32646 <- st_transform(conflict_crs, crs = 32646)
count_32646 <- nrow(conflict_data_32646)
# Count number of conflicts for CRS 32647
conflict_data_32647 <- st_transform(conflict_crs, crs = 32647)
count_32647 <- nrow(conflict_data_32647)

crs_counts <- data.frame(
  CRS = c("EPSG: 32646", "EPSG: 32647"),
  Conflicts_Count = c(count_32646, count_32647)
)

print(crs_counts)
          CRS Conflicts_Count
1 EPSG: 32646           87746
2 EPSG: 32647           87746

Since there is no difference in the count, I will decide to focus on UTM zone 47N (EPSG:32647), east of Myanmar, for the purpose of this exercise. The st_transform() function below converts the CRS of the sf object to EPSG:32647.

# Import boundary data
boundary_sf <- st_read(dsn = "data/geospatial",layer = "mmr_polbnda_adm2_250k_mimu") %>% st_transform(crs = 32647)
Reading layer `mmr_polbnda_adm2_250k_mimu' from data source 
  `C:\SamanthaxFoo\IS415-GAA\Take-home_Ex\Take-home_Ex1\data\geospatial' 
  using driver `ESRI Shapefile'
Simple feature collection with 80 features and 7 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: 92.1721 ymin: 9.696844 xmax: 101.17 ymax: 28.54554
Geodetic CRS:  WGS 84

In the code below, we can notice that the ESPG code has been updated to 32647.

# Check for changes
st_crs(boundary_sf)
Coordinate Reference System:
  User input: EPSG:32647 
  wkt:
PROJCRS["WGS 84 / UTM zone 47N",
    BASEGEOGCRS["WGS 84",
        ENSEMBLE["World Geodetic System 1984 ensemble",
            MEMBER["World Geodetic System 1984 (Transit)"],
            MEMBER["World Geodetic System 1984 (G730)"],
            MEMBER["World Geodetic System 1984 (G873)"],
            MEMBER["World Geodetic System 1984 (G1150)"],
            MEMBER["World Geodetic System 1984 (G1674)"],
            MEMBER["World Geodetic System 1984 (G1762)"],
            MEMBER["World Geodetic System 1984 (G2139)"],
            ELLIPSOID["WGS 84",6378137,298.257223563,
                LENGTHUNIT["metre",1]],
            ENSEMBLEACCURACY[2.0]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433]],
        ID["EPSG",4326]],
    CONVERSION["UTM zone 47N",
        METHOD["Transverse Mercator",
            ID["EPSG",9807]],
        PARAMETER["Latitude of natural origin",0,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8801]],
        PARAMETER["Longitude of natural origin",99,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8802]],
        PARAMETER["Scale factor at natural origin",0.9996,
            SCALEUNIT["unity",1],
            ID["EPSG",8805]],
        PARAMETER["False easting",500000,
            LENGTHUNIT["metre",1],
            ID["EPSG",8806]],
        PARAMETER["False northing",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8807]]],
    CS[Cartesian,2],
        AXIS["(E)",east,
            ORDER[1],
            LENGTHUNIT["metre",1]],
        AXIS["(N)",north,
            ORDER[2],
            LENGTHUNIT["metre",1]],
    USAGE[
        SCOPE["Navigation and medium accuracy spatial referencing."],
        AREA["Between 96°E and 102°E, northern hemisphere between equator and 84°N, onshore and offshore. China. Indonesia. Laos. Malaysia - West Malaysia. Mongolia. Myanmar (Burma). Russian Federation. Thailand."],
        BBOX[0,96,84,102]],
    ID["EPSG",32647]]

Here, I will use the plot() function which plots the geometry of the sf object. The st_geometry() function is used to extract the geometry of the mpsz_sf object which includes the districts of Myanmar as shown below.

par(mar = c(0,0,0,0))
plot(st_geometry(boundary_sf))

3. Data Wrangling

3.1 Fixing Incorrect Datatypes

Recall that the earlier inspection of the conflict_data tibble data frame revealed that the datatype indicated for event date is wrongly labelled as a character instead of a date format.

As such, let’s convert the datatype to the correct ‘date’ format as shown below.

# Convert the datatype for event_date
conflict_data$event_date <- as.Date(conflict_data$event_date, format = "%d %B %Y")

# Check for changes
head(conflict_data)
# A tibble: 6 × 28
  event_id_cnty event_date  year time_precision disorder_type      event_type
  <chr>         <date>     <dbl>          <dbl> <chr>              <chr>     
1 MMR64313      2024-06-30  2024              1 Political violence Battles   
2 MMR64313      2024-06-30  2024              1 Political violence Battles   
3 MMR64320      2024-06-30  2024              1 Political violence Battles   
4 MMR64320      2024-06-30  2024              1 Political violence Battles   
5 MMR64321      2024-06-30  2024              1 Political violence Battles   
6 MMR64321      2024-06-30  2024              1 Political violence Battles   
# ℹ 22 more variables: sub_event_type <chr>, actor1 <chr>, assoc_actor_1 <chr>,
#   inter1 <dbl>, interaction <dbl>, civilian_targeting <chr>, iso <dbl>,
#   region <chr>, country <chr>, admin1 <chr>, admin2 <chr>, admin3 <chr>,
#   location <chr>, latitude <dbl>, longitude <dbl>, geo_precision <dbl>,
#   source <chr>, source_scale <chr>, notes <chr>, fatalities <dbl>,
#   tags <chr>, timestamp <dbl>

3.2 Adding new year_quarter column

We will want to create a new column to indicate the specific year and quarter for each conflict event since the spatial analysis will be done later in a quarterly manner.

Extract year and quarter
conflict_data$year_quarter <- paste0(
  year(conflict_data$event_date), 
  " Q", 
  quarter(conflict_data$event_date)
)

# View the new data column
unique(conflict_data$year_quarter)
 [1] "2024 Q2" "2024 Q1" "2023 Q4" "2023 Q3" "2023 Q2" "2023 Q1" "2022 Q4"
 [8] "2022 Q3" "2022 Q2" "2022 Q1" "2021 Q4" "2021 Q3" "2021 Q2" "2021 Q1"

3.3 Fixing Duplicated Event ID in conflict_data Dataframe

As shown, there are presence of duplicates in our dataframe returned by the duplicated() function.

# Check for duplicates
any(duplicated(conflict_data))
[1] TRUE

Based on the duplicated event ID: MMR64313 for instance. We can observe the two records are of the same political violence event happening between two actors on 30/6/2024, between the People’s Defense Force and Military Forces of Myanmar. Upon further research, these two actors are opposing political parties of Myanmar’s ongoing conflict.

# Inspect an instance of the duplciated event IDs
head(conflict_data,2)
# A tibble: 2 × 29
  event_id_cnty event_date  year time_precision disorder_type      event_type
  <chr>         <date>     <dbl>          <dbl> <chr>              <chr>     
1 MMR64313      2024-06-30  2024              1 Political violence Battles   
2 MMR64313      2024-06-30  2024              1 Political violence Battles   
# ℹ 23 more variables: sub_event_type <chr>, actor1 <chr>, assoc_actor_1 <chr>,
#   inter1 <dbl>, interaction <dbl>, civilian_targeting <chr>, iso <dbl>,
#   region <chr>, country <chr>, admin1 <chr>, admin2 <chr>, admin3 <chr>,
#   location <chr>, latitude <dbl>, longitude <dbl>, geo_precision <dbl>,
#   source <chr>, source_scale <chr>, notes <chr>, fatalities <dbl>,
#   tags <chr>, timestamp <dbl>, year_quarter <chr>
Reflection

Should duplicated data be removed in this analysis?

A single event (e.g. MMR64313) can have duplicated rows with different actor1 values, typically due to counterattacks from opposing sides, leading to different data entries into the conflict_data dataset.

Hence, I will remove duplicated events found in the conflict_data dataframe as long as the rows have the same event ID indicated.

Here, I did another check to ensure there is not more than 2 possible repeated event IDs in the first 20 rows of conflict_data.

Check duplicated events for first 20 rows
duplicate_counts_first_20 <- conflict_data %>%
  slice(1:20) %>%            
  group_by(event_id_cnty) %>% 
  summarize(count = n()) %>%  
  filter(count > 1)         

# View the result
print(duplicate_counts_first_20)
# A tibble: 9 × 2
  event_id_cnty count
  <chr>         <int>
1 MMR64313          2
2 MMR64320          2
3 MMR64321          2
4 MMR64323          2
5 MMR64325          2
6 MMR64326          2
7 MMR64328          2
8 MMR64330          2
9 MMR64331          2

With that checked, I’ll remove the duplicated rows with a repeated Event ID but will retain the actor1 value and append it to a new actor2 column, assuming two actors are involved in every conflict.

Remove duplicated rows
# Load necessary library
library(dplyr)

# Retrieve data of duplicated rows and summarize actor2, assoc_actor_2, and inter2
merged_duplicates <- conflict_data %>%
  filter(duplicated(event_id_cnty) | duplicated(event_id_cnty, fromLast = TRUE)) %>%
  arrange(event_id_cnty) %>%
  group_by(event_id_cnty) %>%
  summarize(
    actor2 = last(actor1)
  )

# Keep rows without duplicates
conflict_data_no_duplicates <- conflict_data %>%
  filter(!duplicated(event_id_cnty))

# Update conflict_data dataframe with new columns from merged_duplicates
conflict_data <- conflict_data_no_duplicates %>%
  left_join(merged_duplicates, by = "event_id_cnty") %>%
  mutate(
    # Ensure actor2 exists after the join, if not fill with actor1
    actor2 = coalesce(actor2, actor1)
  )

# View dataframe
head(conflict_data[c('event_id_cnty','actor1','actor2')])
# A tibble: 6 × 3
  event_id_cnty actor1                                           actor2         
  <chr>         <chr>                                            <chr>          
1 MMR64313      People's Defense Force - Mandalay                Military Force…
2 MMR64320      People's Defense Force - Mandalay                Military Force…
3 MMR64321      People's Defense Force - Mandalay                Military Force…
4 MMR64322      Military Forces of Myanmar (2021-)               Military Force…
5 MMR64323      PKDF MMU: People's Knight Defense Force - Myinmu Military Force…
6 MMR64324      Military Forces of Myanmar (2021-)               Military Force…

We can observe that there are no longer any duplicated event IDs in our conflict_data data frame.

any(duplicated(conflict_data))
[1] FALSE

3.4 Converting Aspatial Data to Simple Feature Format

For the purpose of this exercise, we will want to integrate and analyse aspatial data in a geographic context. I’ll do a check if conflict_data needs to be converted to a sf data frame - if it outputs anything else but sf, then it’s not a simple feature data frame!

class(conflict_data)
[1] "tbl_df"     "tbl"        "data.frame"
Observations

We can see that conflict_data is not a sf data frame. Since a non-simple feature data frame does not have a “geometry” column, we’ll need to convert conflict_data into a simple feature data frame

We can convert conflict_data into a simple feature data frame by using st_as_sf() from the sf package. Addiitionally, we will also need to transform coordinate system from geographic (ESPG: 4326) to projected (ESPG: 32647) using st_transform().

# Convert to simple feature format
conflict_data_sf <- st_as_sf(conflict_data, coords = c("longitude", "latitude"), crs=4326) %>% st_transform(crs = 32647)

# Inspect the changes
glimpse(conflict_data_sf)
Rows: 51,553
Columns: 29
$ event_id_cnty      <chr> "MMR64313", "MMR64320", "MMR64321", "MMR64322", "MM…
$ event_date         <date> 2024-06-30, 2024-06-30, 2024-06-30, 2024-06-30, 20…
$ year               <dbl> 2024, 2024, 2024, 2024, 2024, 2024, 2024, 2024, 202…
$ time_precision     <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, …
$ disorder_type      <chr> "Political violence", "Political violence", "Politi…
$ event_type         <chr> "Battles", "Battles", "Battles", "Strategic develop…
$ sub_event_type     <chr> "Armed clash", "Armed clash", "Armed clash", "Chang…
$ actor1             <chr> "People's Defense Force - Mandalay", "People's Defe…
$ assoc_actor_1      <chr> "MDA - AGF: Madaya - The Authentic Genes Force; SST…
$ inter1             <dbl> 3, 3, 3, 1, 3, 1, 1, 2, 1, 1, 1, 1, 3, 3, 3, 7, 1, …
$ interaction        <dbl> 13, 13, 13, 10, 13, 10, 12, 12, 12, 12, 12, 13, 13,…
$ civilian_targeting <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,…
$ iso                <dbl> 104, 104, 104, 104, 104, 104, 104, 104, 104, 104, 1…
$ region             <chr> "Southeast Asia", "Southeast Asia", "Southeast Asia…
$ country            <chr> "Myanmar", "Myanmar", "Myanmar", "Myanmar", "Myanma…
$ admin1             <chr> "Mandalay", "Mandalay", "Mandalay", "Sagaing", "Sag…
$ admin2             <chr> "Mandalay", "Pyinoolwin", "Pyinoolwin", "Shwebo", "…
$ admin3             <chr> "Patheingyi", "Singu", "Thabeikkyin", "Khin-U", "My…
$ location           <chr> "Aung Tha Pyay", "Pin Lel Gyi", "Thabeikkyin", "Khi…
$ geo_precision      <dbl> 2, 1, 2, 2, 2, 2, 1, 1, 2, 1, 1, 2, 1, 1, 2, 2, 1, …
$ source             <chr> "Democratic Voice of Burma; Irrawaddy", "Irrawaddy"…
$ source_scale       <chr> "National", "National", "National", "Subnational-Na…
$ notes              <chr> "On 30 June 2024, near Aung Tha Pyay village (Pathe…
$ fatalities         <dbl> 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
$ tags               <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,…
$ timestamp          <dbl> 1720552468, 1720552468, 1720552468, 1720552468, 172…
$ year_quarter       <chr> "2024 Q2", "2024 Q2", "2024 Q2", "2024 Q2", "2024 Q…
$ actor2             <chr> "Military Forces of Myanmar (2021-)", "Military For…
$ geometry           <POINT [m]> POINT (214961 2452068), POINT (198303.2 24994…
Observations

Notice that a new column called geometry has been added into the data frame. On the other hand, the longitude and latitude columns have been removed from the data frame.

We can further inspect the newly created ‘geometry’ column of conflict_data_sf

# Retrieve geometry column
st_geometry(conflict_data_sf)
Geometry set for 51553 features 
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: -208804.4 ymin: 1103500 xmax: 640934.5 ymax: 3042960
Projected CRS: WGS 84 / UTM zone 47N
First 5 geometries:
POINT (214961 2452068)
POINT (198303.2 2499463)
POINT (189105.4 2533434)
POINT (160913.9 2522331)
POINT (146213 2428487)
Observations

It consists of 51,533 features consisting of point geometric features where the underlying datum is in WGS 84 format.

To ensure that the coordinate system is correctly updated, we can use the st_crs() function where we observe that the ESPG code is correctly indicated as 32647.

# Check CRS format
st_crs(conflict_data_sf)
Coordinate Reference System:
  User input: EPSG:32647 
  wkt:
PROJCRS["WGS 84 / UTM zone 47N",
    BASEGEOGCRS["WGS 84",
        ENSEMBLE["World Geodetic System 1984 ensemble",
            MEMBER["World Geodetic System 1984 (Transit)"],
            MEMBER["World Geodetic System 1984 (G730)"],
            MEMBER["World Geodetic System 1984 (G873)"],
            MEMBER["World Geodetic System 1984 (G1150)"],
            MEMBER["World Geodetic System 1984 (G1674)"],
            MEMBER["World Geodetic System 1984 (G1762)"],
            MEMBER["World Geodetic System 1984 (G2139)"],
            ELLIPSOID["WGS 84",6378137,298.257223563,
                LENGTHUNIT["metre",1]],
            ENSEMBLEACCURACY[2.0]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433]],
        ID["EPSG",4326]],
    CONVERSION["UTM zone 47N",
        METHOD["Transverse Mercator",
            ID["EPSG",9807]],
        PARAMETER["Latitude of natural origin",0,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8801]],
        PARAMETER["Longitude of natural origin",99,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8802]],
        PARAMETER["Scale factor at natural origin",0.9996,
            SCALEUNIT["unity",1],
            ID["EPSG",8805]],
        PARAMETER["False easting",500000,
            LENGTHUNIT["metre",1],
            ID["EPSG",8806]],
        PARAMETER["False northing",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8807]]],
    CS[Cartesian,2],
        AXIS["(E)",east,
            ORDER[1],
            LENGTHUNIT["metre",1]],
        AXIS["(N)",north,
            ORDER[2],
            LENGTHUNIT["metre",1]],
    USAGE[
        SCOPE["Navigation and medium accuracy spatial referencing."],
        AREA["Between 96°E and 102°E, northern hemisphere between equator and 84°N, onshore and offshore. China. Indonesia. Laos. Malaysia - West Malaysia. Mongolia. Myanmar (Burma). Russian Federation. Thailand."],
        BBOX[0,96,84,102]],
    ID["EPSG",32647]]

3.5 Reduce Data File Size

In this section, I will reduce the current Myanmar armed conflict dataset as the time taken for computing the kernel density estimates can take up to 30 minutes long which is not computationally efficient.

1) Remove ‘Protests’ and ‘Riots’ Event Types

I will remove rows in the conflicts_data_sf dataset that don’t focus on the four main event types (Battles, Explosion/Remote violence, Strategic developments, and Violence against civilians), as mentioned in the exercise brief.

conflict_data_sf <- conflict_data_sf %>%
  filter(!(event_type %in% c("Protests", "Riots")))

unique(conflict_data_sf$event_type)
[1] "Battles"                    "Strategic developments"    
[3] "Violence against civilians" "Explosions/Remote violence"

2) Remove unused columns in boundary_sf

As seen, there are 8 columns in the simple feature data frame of boundary_sf.

# Inspect first rows of data in boundary_sf
head(boundary_sf)
Simple feature collection with 6 features and 7 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -14915.04 ymin: 1736124 xmax: 187961.7 ymax: 2051144
Projected CRS: WGS 84 / UTM zone 47N
  OBJECTID         ST ST_PCODE        DT   DT_PCODE      DT_MMR PCode_V
1        1 Ayeyarwady   MMR017  Hinthada MMR017D002    ဟင်္သာတခရိုင်     9.4
2        2 Ayeyarwady   MMR017   Labutta MMR017D004    လပွတ္တာခရိုင်     9.4
3        3 Ayeyarwady   MMR017    Maubin MMR017D005     မအူပင်ခရိုင်     9.4
4        4 Ayeyarwady   MMR017 Myaungmya MMR017D003 မြောင်းမြခရိုင်     9.4
5        5 Ayeyarwady   MMR017   Pathein MMR017D001      ပုသိမ်ခရိုင်     9.4
6        6 Ayeyarwady   MMR017    Pyapon MMR017D006     ဖျာပုံခရိုင်     9.4
                        geometry
1 MULTIPOLYGON (((90859.89 20...
2 MULTIPOLYGON (((75991.51 17...
3 MULTIPOLYGON (((115559 1928...
4 MULTIPOLYGON (((39919.39 18...
5 MULTIPOLYGON (((-6302.348 1...
6 MULTIPOLYGON (((93411.72 17...

I will remove ’DT_MMR” column as we already have the District Name in English in DT and won’t require the district names in Myanmar Language. Next, we will remove the coded versions of ST (state/region) and DT (district) columns, namely ST_PCODE and DT_PCODE. Additionally, we won’t need the PCode_V column since we will be dropping the PCODE column too.

boundary_sf <- boundary_sf %>% dplyr::select('OBJECTID', 'ST', 'DT','geometry')
summary(boundary_sf)
    OBJECTID          ST                 DT                     geometry 
 Min.   : 1.00   Length:80          Length:80          MULTIPOLYGON :80  
 1st Qu.:20.75   Class :character   Class :character   epsg:32647   : 0  
 Median :40.50   Mode  :character   Mode  :character   +proj=utm ...: 0  
 Mean   :40.50                                                           
 3rd Qu.:60.25                                                           
 Max.   :80.00                                                           

3) Remove unused columns in conflict_data

I will also remove unnecessary columns of the conflict_data data frame that won’t be used in our spatial analysis later. I’ll rename admin1 to ST (states) and admin2 to DT (districts) for easier reference.

Remove unnecessary columns
conflict_data_sf <- conflict_data_sf %>%
  select(event_id_cnty, event_date, year_quarter, disorder_type, 
         event_type, location, admin1, admin2, geometry, actor1, actor2,
         interaction, fatalities) %>%
  rename(
    ST = admin1,
    DT = admin2
  )

summary(conflict_data_sf)
 event_id_cnty        event_date         year_quarter       disorder_type     
 Length:42608       Min.   :2021-01-01   Length:42608       Length:42608      
 Class :character   1st Qu.:2022-01-10   Class :character   Class :character  
 Mode  :character   Median :2022-10-13   Mode  :character   Mode  :character  
                    Mean   :2022-10-29                                        
                    3rd Qu.:2023-08-29                                        
                    Max.   :2024-06-30                                        
  event_type          location              ST                 DT           
 Length:42608       Length:42608       Length:42608       Length:42608      
 Class :character   Class :character   Class :character   Class :character  
 Mode  :character   Mode  :character   Mode  :character   Mode  :character  
                                                                            
                                                                            
                                                                            
          geometry        actor1             actor2           interaction   
 POINT        :42608   Length:42608       Length:42608       Min.   :10.00  
 epsg:32647   :    0   Class :character   Class :character   1st Qu.:13.00  
 +proj=utm ...:    0   Mode  :character   Mode  :character   Median :17.00  
                                                             Mean   :18.86  
                                                             3rd Qu.:17.00  
                                                             Max.   :80.00  
   fatalities    
 Min.   :  0.00  
 1st Qu.:  0.00  
 Median :  0.00  
 Mean   :  1.27  
 3rd Qu.:  1.00  
 Max.   :201.00  

3.6 Converting Simple Features Data Frame into ppp Object

It is important that we convert conflict_data_sf (a simple feature data frame) into a planer point pattern (ppp) object format, since the spatstat package that we’ll be using for the Spatial Point Pattern Analysis later is specifically designed for working with ppp-formated data. Additionally, I will begin with categorising the ppp objects into their unique year_quarter category.

Create ppp objects based on year_quarter category
# Create an empty list to store the ppp objects
ppp_list <- list()

# Loop through each unique year_quarter category
for (yq in unique(conflict_data_sf$year_quarter)) {
  # Subset the data for the current year_quarter
  subset_data_sf <- conflict_data_sf %>% filter(year_quarter == yq)
  
  # Convert the subset to a ppp object
  subset_ppp <- as.ppp(subset_data_sf$geometry)
  
  # Add the ppp object to the list
  ppp_list[[yq]] <- subset_ppp
}

# Check list
ppp_list
$`2024 Q2`
Planar point pattern: 2788 points
window: rectangle = [-208804.4, 597543.7] x [1103500.1, 3026504.9] units

$`2024 Q1`
Planar point pattern: 3186 points
window: rectangle = [-207135, 591875.9] x [1245380, 3026504.9] units

$`2023 Q4`
Planar point pattern: 3627 points
window: rectangle = [-206931.7, 604775.1] x [1103500.1, 3020772.2] units

$`2023 Q3`
Planar point pattern: 3010 points
window: rectangle = [-197883.4, 518300.4] x [1103500.1, 3027041.8] units

$`2023 Q2`
Planar point pattern: 2745 points
window: rectangle = [-191261.5, 518300.4] x [1103500.1, 3006372.9] units

$`2023 Q1`
Planar point pattern: 3101 points
window: rectangle = [-199243.8, 591875.9] x [1103500.1, 3026504.9] units

$`2022 Q4`
Planar point pattern: 3296 points
window: rectangle = [-206531.5, 518300.4] x [1103500.1, 2931517.1] units

$`2022 Q3`
Planar point pattern: 3486 points
window: rectangle = [-206196.6, 568361.5] x [1103500.1, 3026504.9] units

$`2022 Q2`
Planar point pattern: 3580 points
window: rectangle = [-206931.7, 640934.5] x [1103500.1, 3026504.9] units

$`2022 Q1`
Planar point pattern: 3563 points
window: rectangle = [-204784, 591875.9] x [1103500.1, 3026504.9] units

$`2021 Q4`
Planar point pattern: 3844 points
window: rectangle = [-200024.3, 591875.9] x [1103500.1, 3042960.3] units

$`2021 Q3`
Planar point pattern: 2754 points
window: rectangle = [-193181.1, 591875.9] x [1103500.1, 3042960.3] units

$`2021 Q2`
Planar point pattern: 2916 points
window: rectangle = [-191409.1, 640934.5] x [1132472.1, 3042960.3] units

$`2021 Q1`
Planar point pattern: 712 points
window: rectangle = [-203795.3, 591875.9] x [1375186.1, 3026504.9] units

We can visualise the spread of conflict events across each quarter from January 2021 to June 2024 using the plot() function as shown below.

Visualise the spread of conflicts by year_quarter
# Ensure 'year_quarter' is a factor
conflict_data_sf$year_quarter <- as.factor(conflict_data_sf$year_quarter)

# Loop through each unique year_quarter and create separate plots
year_quarters <- unique(conflict_data_sf$year_quarter)

# Set up a grid layout for multiple plots (adjust 'mfrow' as needed)
par(mfrow = c(4,4))
par(mar = c(0,0,1,0))

# Loop through each year_quarter and plot
for (yq in year_quarters) {
  subset_data_sf <- conflict_data_sf[conflict_data_sf$year_quarter == yq, ]
  conflict_data_ppp <- as.ppp(subset_data_sf$geometry)
  
  # Plot each subset ppp object
  plot(conflict_data_ppp, main = paste(yq))
}

Observations

It is noticeable that there conflict events have occured more frequently since 2021 as points plotted on the graph have gotten darker across 2021 to 2024. We can also observe the possibility of duplicated events occurring from the darker spots in the plot, in which it appears more intense in Myanmar’s central and west regions.

3.7 Creating owin object

When analysing spatial point patterns, it is a good practice to confine the analysis with a geographical area, that is Myanmar’s boundary in this case. In spatstat, an object called owin is specially designed to represent this polygonal region.

The code chunk below is used to convert the boundary_data_sf simple feature data frame into an owin object of spatstat.

# Convert to owin object
myanmar_owin <- as.owin(boundary_sf)

# Visualise the owin object
plot(myanmar_owin)

Observations

From my observations, the as.owin() function converts the boundary_data_sf spatial boundary into a window object that represents the outer boundary of the spatial region and does not handle internal structures or districts we previously saw from the plot of boundary_data_sf.

We can also take a quick look at the owin object properties as shown. I will be converting it to a data frame for the purposes of getting a quick glimpse of the object.

# Summary info of owin object
owin_df <- as.data.frame(myanmar_owin)
print(head(owin_df))
         x       y id sign
1 56519.39 2741919  1   -1
2 56917.28 2741947  1   -1
3 57000.15 2741973  1   -1
4 57068.51 2741994  1   -1
5 57221.44 2742142  1   -1
6 57068.51 2741994  1   -1

3.8 Combining ppp Object and owin Object

In this last step of geospatial data wrangling, I will mask all ppp object with the owin object I created earlier to put in place all conflict events within the boundary of Myanmar. Doing so can also optimise the memory usage for large datasets. I’ll also make the density values more comprehensible by re-scaling the density values from metres to kilometres using rescale().

# Initialize an empty list to store masked ppp objects
masked_ppp_list_km <- list()

# Iterate over each ppp object in the list
for (quarter in names(ppp_list)) {
  ppp_obj <- ppp_list[[quarter]]
  # Mask the ppp object with the owin object
  masked_ppp <- ppp_obj[myanmar_owin]
  masked_ppp_km <- rescale(masked_ppp, 1000, "km")
  # Store the masked ppp object in the new list
  masked_ppp_list_km[[quarter]] <- masked_ppp_km
}
# Inspect 2024 Q2 masked ppp object
summary(masked_ppp_list_km$`2024 Q2`)
Planar point pattern:  2788 points
Average intensity 0.004162974 points per square km

*Pattern contains duplicated points*

Coordinates are given to 16 decimal places

Window: polygonal boundary
1077 separate polygons (510 holes)
                    vertices         area relative.area
polygon 1 (hole)           9 -4.40155e-08     -6.57e-14
polygon 2 (hole)           3 -4.82818e-08     -7.21e-14
polygon 3 (hole)           3 -2.78345e-09     -4.16e-15
polygon 4 (hole)           3 -4.99415e-09     -7.46e-15
polygon 5 (hole)           3 -6.70916e-10     -1.00e-15
polygon 6 (hole)           3 -4.99758e-09     -7.46e-15
polygon 7 (hole)           6 -3.24429e-09     -4.84e-15
polygon 8 (hole)          11 -5.35507e-08     -8.00e-14
polygon 9 (hole)           3 -3.06875e-09     -4.58e-15
polygon 10 (hole)          3 -7.63672e-09     -1.14e-14
polygon 11 (hole)          3 -8.61828e-11     -1.29e-16
polygon 12 (hole)          3 -2.44104e-08     -3.64e-14
polygon 13 (hole)          3 -4.88141e-11     -7.29e-17
polygon 14 (hole)          4 -1.74466e-10     -2.61e-16
polygon 15 (hole)          4 -2.30890e-11     -3.45e-17
polygon 16 (hole)          4 -6.25346e-14     -9.34e-20
polygon 17 (hole)          4 -3.68040e-08     -5.50e-14
polygon 18 (hole)          3 -1.27228e-12     -1.90e-18
polygon 19 (hole)          4 -1.25504e-09     -1.87e-15
polygon 20 (hole)          3 -1.56980e-08     -2.34e-14
polygon 21 (hole)          4 -1.13609e-12     -1.70e-18
polygon 22 (hole)          4 -8.02868e-08     -1.20e-13
polygon 23 (hole)          4 -3.98057e-11     -5.94e-17
polygon 24 (hole)          4 -2.89739e-08     -4.33e-14
polygon 25 (hole)          5 -7.39370e-08     -1.10e-13
polygon 26 (hole)          6 -1.08246e-07     -1.62e-13
polygon 27 (hole)         13 -3.27066e-07     -4.88e-13
polygon 28 (hole)          3 -3.55050e-09     -5.30e-15
polygon 29 (hole)          3 -2.84064e-08     -4.24e-14
polygon 30 (hole)          3 -8.07641e-10     -1.21e-15
polygon 31 (hole)          4 -5.20948e-09     -7.78e-15
polygon 32 (hole)          4 -3.21399e-10     -4.80e-16
polygon 33 (hole)          4 -2.54955e-11     -3.81e-17
polygon 34 (hole)          3 -4.29638e-10     -6.42e-16
polygon 35 (hole)          6 -6.02383e-10     -8.99e-16
polygon 36 (hole)          4 -4.91852e-10     -7.34e-16
polygon 37                26  2.85778e+00      4.27e-06
polygon 38 (hole)          3 -6.80850e-13     -1.02e-18
polygon 39 (hole)          4 -3.09085e-13     -4.62e-19
polygon 40 (hole)          6 -3.61049e-12     -5.39e-18
polygon 41 (hole)          7 -5.33310e-12     -7.96e-18
polygon 42 (hole)          4 -3.62663e-12     -5.42e-18
polygon 43 (hole)          3 -2.00302e-14     -2.99e-20
polygon 44 (hole)          3 -2.85255e-12     -4.26e-18
polygon 45 (hole)          4 -1.15535e-12     -1.73e-18
polygon 46 (hole)         17 -6.25657e-13     -9.34e-19
polygon 47 (hole)         20 -1.53263e-11     -2.29e-17
polygon 48 (hole)          3 -1.42406e-14     -2.13e-20
polygon 49                 3  0.00000e+00      0.00e+00
polygon 50 (hole)          4 -2.46048e-12     -3.67e-18
polygon 51 (hole)          4 -2.93672e-13     -4.39e-19
polygon 52 (hole)          3 -1.35167e-16     -2.02e-22
polygon 53 (hole)          3 -3.66243e-13     -5.47e-19
polygon 54 (hole)         13 -3.22774e-12     -4.82e-18
polygon 55 (hole)          7 -3.79713e-12     -5.67e-18
polygon 56 (hole)          4 -6.06083e-13     -9.05e-19
polygon 57 (hole)          6 -5.94835e-12     -8.88e-18
polygon 58 (hole)          4 -4.11305e-13     -6.14e-19
polygon 59 (hole)         12 -7.16981e-12     -1.07e-17
polygon 60 (hole)          4 -1.96241e-14     -2.93e-20
polygon 61 (hole)          4 -6.16304e-12     -9.20e-18
polygon 62 (hole)          4 -2.13366e-13     -3.19e-19
polygon 63 (hole)          4 -2.21852e-12     -3.31e-18
polygon 64 (hole)          3 -4.22166e-15     -6.30e-21
polygon 65 (hole)          5 -6.21282e-13     -9.28e-19
polygon 66 (hole)          4 -5.02821e-16     -7.51e-22
polygon 67 (hole)          4 -3.30895e-13     -4.94e-19
polygon 68 (hole)         17 -1.19462e-11     -1.78e-17
polygon 69 (hole)          3 -1.27084e-13     -1.90e-19
polygon 70 (hole)          4 -1.80406e-13     -2.69e-19
polygon 71 (hole)          8 -3.87178e-12     -5.78e-18
polygon 72 (hole)          4 -7.14736e-13     -1.07e-18
polygon 73 (hole)          5 -3.43923e-12     -5.14e-18
polygon 74                 3  0.00000e+00      0.00e+00
polygon 75 (hole)          8 -2.57732e-12     -3.85e-18
polygon 76 (hole)          4 -1.36629e-12     -2.04e-18
polygon 77 (hole)         13 -1.46870e-11     -2.19e-17
polygon 78 (hole)          3 -3.08778e-14     -4.61e-20
polygon 79 (hole)         18 -1.78311e-12     -2.66e-18
polygon 80 (hole)         28 -1.67874e-11     -2.51e-17
polygon 81 (hole)          4 -4.46113e-13     -6.66e-19
polygon 82 (hole)          3 -8.98160e-15     -1.34e-20
polygon 83 (hole)          4 -1.43028e-12     -2.14e-18
polygon 84 (hole)          4 -6.38438e-14     -9.53e-20
polygon 85 (hole)         19 -1.89999e-11     -2.84e-17
polygon 86 (hole)          4 -1.49124e-13     -2.23e-19
polygon 87 (hole)          3 -3.68533e-14     -5.50e-20
polygon 88 (hole)          4 -4.45374e-13     -6.65e-19
polygon 89 (hole)          4 -4.83479e-13     -7.22e-19
polygon 90 (hole)         12 -9.32407e-12     -1.39e-17
polygon 91 (hole)          3 -8.84897e-13     -1.32e-18
polygon 92 (hole)         12 -8.41711e-12     -1.26e-17
polygon 93 (hole)          3 -3.73405e-14     -5.58e-20
polygon 94 (hole)          3 -1.40773e-12     -2.10e-18
polygon 95 (hole)          4 -7.84182e-13     -1.17e-18
polygon 96 (hole)          5 -4.10716e-13     -6.13e-19
polygon 97 (hole)         10 -5.89467e-12     -8.80e-18
polygon 98 (hole)          3 -3.54370e-15     -5.29e-21
polygon 99 (hole)         11 -6.90903e-12     -1.03e-17
polygon 100 (hole)         3 -1.89866e-13     -2.84e-19
polygon 101 (hole)         3 -3.09480e-15     -4.62e-21
polygon 102 (hole)         3 -3.95637e-17     -5.91e-23
polygon 103 (hole)         4 -1.88731e-16     -2.82e-22
polygon 104 (hole)         3 -5.93215e-17     -8.86e-23
polygon 105 (hole)         6 -1.10389e-12     -1.65e-18
polygon 106 (hole)         5 -1.62107e-13     -2.42e-19
polygon 107 (hole)         4 -6.07269e-13     -9.07e-19
polygon 108 (hole)         3 -4.17062e-14     -6.23e-20
polygon 109 (hole)         7 -4.10998e-12     -6.14e-18
polygon 110 (hole)         7 -2.99640e-12     -4.47e-18
polygon 111 (hole)         4 -8.08268e-13     -1.21e-18
polygon 112 (hole)         8 -1.72022e-12     -2.57e-18
polygon 113 (hole)         3 -5.46853e-13     -8.17e-19
polygon 114 (hole)         3 -9.21233e-15     -1.38e-20
polygon 115 (hole)         3 -2.34233e-14     -3.50e-20
polygon 116 (hole)         9 -1.27183e-12     -1.90e-18
polygon 117 (hole)         4 -8.58065e-13     -1.28e-18
polygon 118 (hole)         6 -4.85387e-12     -7.25e-18
polygon 119 (hole)         6 -4.38024e-12     -6.54e-18
polygon 120 (hole)         6 -1.18459e-12     -1.77e-18
polygon 121 (hole)        18 -5.24919e-12     -7.84e-18
polygon 122                3  0.00000e+00      0.00e+00
polygon 123 (hole)        10 -3.29610e-12     -4.92e-18
polygon 124 (hole)         4 -8.31786e-15     -1.24e-20
polygon 125 (hole)         8 -5.92159e-12     -8.84e-18
polygon 126 (hole)         4 -2.92352e-15     -4.37e-21
polygon 127 (hole)         3 -2.22606e-13     -3.32e-19
polygon 128 (hole)        16 -1.60557e-11     -2.40e-17
polygon 129 (hole)        10 -4.89387e-12     -7.31e-18
polygon 130 (hole)         4 -1.48133e-13     -2.21e-19
polygon 131 (hole)         4 -2.04528e-12     -3.05e-18
polygon 132 (hole)         4 -4.07107e-13     -6.08e-19
polygon 133 (hole)         4 -2.13923e-12     -3.19e-18
polygon 134 (hole)         4 -2.58593e-17     -3.86e-23
polygon 135 (hole)         3 -8.40799e-14     -1.26e-19
polygon 136 (hole)         4 -5.03212e-13     -7.51e-19
polygon 137 (hole)         4 -1.10455e-12     -1.65e-18
polygon 138 (hole)         4 -3.74743e-14     -5.60e-20
polygon 139 (hole)         6 -1.42384e-12     -2.13e-18
polygon 140 (hole)         4 -1.96713e-14     -2.94e-20
polygon 141 (hole)         4 -7.29770e-13     -1.09e-18
polygon 142 (hole)         3 -7.99677e-16     -1.19e-21
polygon 143 (hole)         3 -1.00370e-15     -1.50e-21
polygon 144 (hole)         3 -9.07286e-16     -1.35e-21
polygon 145 (hole)         4 -2.63668e-14     -3.94e-20
polygon 146 (hole)         5 -8.51993e-13     -1.27e-18
polygon 147 (hole)         4 -3.15631e-13     -4.71e-19
polygon 148 (hole)         4 -7.89181e-13     -1.18e-18
polygon 149 (hole)         6 -4.24179e-13     -6.33e-19
polygon 150 (hole)        16 -1.00987e-11     -1.51e-17
polygon 151 (hole)         3 -6.53241e-18     -9.75e-24
polygon 152 (hole)         6 -2.87035e-12     -4.29e-18
polygon 153 (hole)         3 -2.64530e-15     -3.95e-21
polygon 154 (hole)         8 -1.99962e-12     -2.99e-18
polygon 155 (hole)         3 -2.45345e-13     -3.66e-19
polygon 156 (hole)         3 -3.99767e-14     -5.97e-20
polygon 157 (hole)         6 -1.08025e-12     -1.61e-18
polygon 158 (hole)         3 -1.22938e-13     -1.84e-19
polygon 159 (hole)         5 -3.07870e-13     -4.60e-19
polygon 160 (hole)        12 -4.81197e-12     -7.19e-18
polygon 161 (hole)         4 -8.16819e-14     -1.22e-19
polygon 162 (hole)         3 -1.28714e-12     -1.92e-18
polygon 163 (hole)         4 -1.87938e-13     -2.81e-19
polygon 164 (hole)        11 -1.79968e-12     -2.69e-18
polygon 165 (hole)         3 -5.93827e-15     -8.87e-21
polygon 166 (hole)         3 -4.50189e-13     -6.72e-19
polygon 167 (hole)         4 -4.87377e-13     -7.28e-19
polygon 168 (hole)        13 -1.01705e-11     -1.52e-17
polygon 169 (hole)         3 -7.37542e-13     -1.10e-18
polygon 170 (hole)         4 -3.62179e-13     -5.41e-19
polygon 171                3  9.13135e-19      1.36e-24
polygon 172 (hole)         8 -1.74082e-12     -2.60e-18
polygon 173 (hole)         3 -9.31514e-14     -1.39e-19
polygon 174               43  7.32477e+00      1.09e-05
polygon 175 (hole)         3 -6.79111e-14     -1.01e-19
polygon 176 (hole)        19 -1.27551e-11     -1.90e-17
polygon 177 (hole)         3 -7.65776e-15     -1.14e-20
polygon 178 (hole)         9 -4.34371e-12     -6.49e-18
polygon 179 (hole)         4 -7.09328e-15     -1.06e-20
polygon 180 (hole)         6 -8.43021e-12     -1.26e-17
polygon 181 (hole)        30 -5.76461e-12     -8.61e-18
polygon 182 (hole)        12 -5.45083e-12     -8.14e-18
polygon 183 (hole)        16 -1.90750e-11     -2.85e-17
polygon 184 (hole)        15 -2.10315e-11     -3.14e-17
polygon 185 (hole)        52 -3.02970e-11     -4.52e-17
polygon 186              103  1.86991e+01      2.79e-05
polygon 187 (hole)        10 -7.40529e-12     -1.11e-17
polygon 188 (hole)         6 -2.03062e-12     -3.03e-18
polygon 189 (hole)         7 -1.02421e-12     -1.53e-18
polygon 190 (hole)         3 -1.32600e-19     -1.98e-25
polygon 191 (hole)         3 -5.60223e-13     -8.37e-19
polygon 192 (hole)        10 -2.37610e-12     -3.55e-18
polygon 193 (hole)         3 -1.90831e-16     -2.85e-22
polygon 194 (hole)         4 -1.76514e-12     -2.64e-18
polygon 195 (hole)         3 -1.37598e-16     -2.05e-22
polygon 196 (hole)         4 -5.28260e-13     -7.89e-19
polygon 197               37  9.32316e+00      1.39e-05
polygon 198 (hole)         4 -2.33117e-13     -3.48e-19
polygon 199 (hole)         4 -6.59209e-13     -9.84e-19
polygon 200 (hole)         4 -5.12582e-13     -7.65e-19
polygon 201 (hole)         4 -4.70065e-13     -7.02e-19
polygon 202 (hole)         3 -2.10303e-13     -3.14e-19
polygon 203 (hole)         6 -1.36386e-12     -2.04e-18
polygon 204 (hole)         4 -2.18244e-13     -3.26e-19
polygon 205 (hole)         4 -3.72475e-13     -5.56e-19
polygon 206 (hole)         4 -8.32022e-13     -1.24e-18
polygon 207 (hole)         4 -1.08153e-12     -1.61e-18
polygon 208 (hole)         3 -5.92333e-17     -8.84e-23
polygon 209 (hole)         7 -6.03242e-12     -9.01e-18
polygon 210 (hole)         6 -2.27756e-12     -3.40e-18
polygon 211 (hole)         4 -1.53368e-12     -2.29e-18
polygon 212 (hole)         3 -2.33223e-14     -3.48e-20
polygon 213 (hole)         6 -8.61964e-12     -1.29e-17
polygon 214              371  2.43869e+02      3.64e-04
polygon 215 (hole)         4 -8.45587e-13     -1.26e-18
polygon 216              297  2.84905e+02      4.25e-04
polygon 217 (hole)         4 -3.40250e-13     -5.08e-19
polygon 218 (hole)         6 -1.39288e-12     -2.08e-18
polygon 219 (hole)         4 -1.27223e-13     -1.90e-19
polygon 220 (hole)         4 -9.59107e-13     -1.43e-18
polygon 221 (hole)         8 -1.39768e-11     -2.09e-17
polygon 222 (hole)         3 -3.34368e-14     -4.99e-20
polygon 223 (hole)         4 -9.58519e-13     -1.43e-18
polygon 224 (hole)         5 -1.32345e-13     -1.98e-19
polygon 225 (hole)         4 -5.98296e-13     -8.93e-19
polygon 226 (hole)         7 -3.05539e-12     -4.56e-18
polygon 227 (hole)         3 -1.80098e-15     -2.69e-21
polygon 228 (hole)         9 -6.38785e-12     -9.54e-18
polygon 229               33  1.68222e+01      2.51e-05
polygon 230 (hole)         8 -1.04272e-11     -1.56e-17
polygon 231               33  4.47665e-01      6.68e-07
polygon 232               19  1.34593e-01      2.01e-07
polygon 233               39  1.36327e+00      2.04e-06
polygon 234              137  1.55547e+02      2.32e-04
polygon 235               36  8.76479e+00      1.31e-05
polygon 236               79  3.08116e+01      4.60e-05
polygon 237              388  2.25271e+02      3.36e-04
polygon 238              316  7.78512e+01      1.16e-04
polygon 239               13  1.09564e-01      1.64e-07
polygon 240               18  3.49727e-01      5.22e-07
polygon 241               31  1.23017e+00      1.84e-06
polygon 242               16  6.55537e-01      9.79e-07
polygon 243               24  8.49487e-01      1.27e-06
polygon 244               30  2.54436e+00      3.80e-06
polygon 245              336  4.15806e+01      6.21e-05
polygon 246              330  1.69190e+02      2.53e-04
polygon 247               47  1.08035e+01      1.61e-05
polygon 248               39  4.94369e+00      7.38e-06
polygon 249               23  2.72438e+00      4.07e-06
polygon 250               33  5.70263e+00      8.52e-06
polygon 251               90  4.20329e+01      6.28e-05
polygon 252               28  1.35341e+00      2.02e-06
polygon 253              225  1.08816e+02      1.62e-04
polygon 254               33  9.16670e+00      1.37e-05
polygon 255              192  7.02655e+01      1.05e-04
polygon 256               49  1.49245e+01      2.23e-05
polygon 257               98  1.79076e+01      2.67e-05
polygon 258                6  6.37552e-01      9.52e-07
polygon 259               49  1.01233e+01      1.51e-05
polygon 260              141  3.43053e+01      5.12e-05
polygon 261              195  3.24345e+01      4.84e-05
polygon 262               51  3.38313e+00      5.05e-06
polygon 263               34  2.01400e+00      3.01e-06
polygon 264               13  2.50435e-01      3.74e-07
polygon 265                9  9.04824e-02      1.35e-07
polygon 266               34  4.61794e+00      6.90e-06
polygon 267               17  4.58200e-01      6.84e-07
polygon 268               15  2.74776e-01      4.10e-07
polygon 269               21  5.34978e-01      7.99e-07
polygon 270               19  4.55347e-01      6.80e-07
polygon 271               71  3.42557e+00      5.11e-06
polygon 272               24  1.32420e+00      1.98e-06
polygon 273               15  3.26247e-01      4.87e-07
polygon 274               39  8.65790e-01      1.29e-06
polygon 275               43  1.41627e+00      2.11e-06
polygon 276               24  7.52068e-01      1.12e-06
polygon 277               96  1.32101e+01      1.97e-05
polygon 278               38  1.18003e+00      1.76e-06
polygon 279              429  5.99087e+02      8.95e-04
polygon 280               13  1.74105e-01      2.60e-07
polygon 281               19  2.52336e-01      3.77e-07
polygon 282               16  3.11495e-01      4.65e-07
polygon 283               11  9.11047e-02      1.36e-07
polygon 284               12  2.13470e-01      3.19e-07
polygon 285               17  5.82663e-01      8.70e-07
polygon 286               56  2.60440e+01      3.89e-05
polygon 287              107  4.91389e+00      7.34e-06
polygon 288               51  2.79076e+00      4.17e-06
polygon 289               89  1.61156e+01      2.41e-05
polygon 290               28  1.30499e+00      1.95e-06
polygon 291               11  1.27616e-01      1.91e-07
polygon 292               34  2.54199e+00      3.80e-06
polygon 293               27  1.72476e+00      2.58e-06
polygon 294               37  2.01882e+00      3.01e-06
polygon 295               23  1.65571e+00      2.47e-06
polygon 296               33  3.05816e+00      4.57e-06
polygon 297               14  3.23153e-01      4.83e-07
polygon 298               91  1.51209e+01      2.26e-05
polygon 299               12  2.42901e-01      3.63e-07
polygon 300               11  1.37889e-01      2.06e-07
polygon 301               58  2.29751e+01      3.43e-05
polygon 302               48  5.10265e+00      7.62e-06
polygon 303               22  1.30706e+00      1.95e-06
polygon 304               15  3.49480e-01      5.22e-07
polygon 305               17  1.57570e+00      2.35e-06
polygon 306               34  3.68725e+00      5.51e-06
polygon 307               34  5.21904e+00      7.79e-06
polygon 308               24  5.42734e+00      8.10e-06
polygon 309              422  4.66497e+02      6.97e-04
polygon 310              142  2.98767e+01      4.46e-05
polygon 311              132  2.18707e+01      3.27e-05
polygon 312               19  5.88230e-01      8.78e-07
polygon 313               22  1.77611e+00      2.65e-06
polygon 314 (hole)         4 -6.42926e-10     -9.60e-16
polygon 315               40  4.09952e+00      6.12e-06
polygon 316 (hole)         3 -2.52076e-09     -3.76e-15
polygon 317 (hole)        24 -5.61620e-07     -8.39e-13
polygon 318               28  1.47685e+00      2.21e-06
polygon 319 (hole)         8 -8.93461e-08     -1.33e-13
polygon 320 (hole)         8 -1.34123e-07     -2.00e-13
polygon 321               67  9.99685e+00      1.49e-05
polygon 322 (hole)         8 -6.36345e-08     -9.50e-14
polygon 323 (hole)         6 -7.57978e-08     -1.13e-13
polygon 324 (hole)         3 -4.32022e-12     -6.45e-18
polygon 325 (hole)         4 -2.24224e-08     -3.35e-14
polygon 326 (hole)         5 -9.35663e-08     -1.40e-13
polygon 327 (hole)        13 -1.39251e-07     -2.08e-13
polygon 328 (hole)         7 -1.46533e-07     -2.19e-13
polygon 329 (hole)         6 -7.59071e-08     -1.13e-13
polygon 330 (hole)         7 -5.15915e-09     -7.70e-15
polygon 331 (hole)        10 -4.09465e-08     -6.11e-14
polygon 332 (hole)        21 -1.86624e-07     -2.79e-13
polygon 333 (hole)         3 -1.69684e-08     -2.53e-14
polygon 334 (hole)         8 -4.08363e-08     -6.10e-14
polygon 335 (hole)         3 -1.73631e-08     -2.59e-14
polygon 336 (hole)         5 -1.41820e-08     -2.12e-14
polygon 337 (hole)        10 -6.15226e-08     -9.19e-14
polygon 338 (hole)         8 -2.90782e-08     -4.34e-14
polygon 339 (hole)         8 -6.15579e-08     -9.19e-14
polygon 340 (hole)         4 -7.46371e-10     -1.11e-15
polygon 341 (hole)         4 -1.05432e-08     -1.57e-14
polygon 342 (hole)         4 -8.25148e-09     -1.23e-14
polygon 343 (hole)         3 -2.59991e-08     -3.88e-14
polygon 344 (hole)         3 -9.09758e-09     -1.36e-14
polygon 345 (hole)         3 -1.96869e-08     -2.94e-14
polygon 346 (hole)         3 -4.01481e-11     -5.99e-17
polygon 347 (hole)         5 -1.49421e-08     -2.23e-14
polygon 348 (hole)         4 -9.23861e-09     -1.38e-14
polygon 349                3  4.91147e-16      7.33e-22
polygon 350 (hole)         7 -5.73304e-08     -8.56e-14
polygon 351 (hole)         4 -4.26735e-08     -6.37e-14
polygon 352 (hole)         7 -1.69748e-08     -2.53e-14
polygon 353 (hole)         3 -1.05617e-08     -1.58e-14
polygon 354               25  4.82266e-01      7.20e-07
polygon 355 (hole)         6 -3.72405e-08     -5.56e-14
polygon 356 (hole)         7 -3.08332e-08     -4.60e-14
polygon 357 (hole)         4 -6.69857e-09     -1.00e-14
polygon 358 (hole)         4 -1.15937e-08     -1.73e-14
polygon 359 (hole)         4 -3.11046e-09     -4.64e-15
polygon 360 (hole)        12 -1.45689e-07     -2.18e-13
polygon 361 (hole)         4 -5.64800e-08     -8.43e-14
polygon 362 (hole)         3 -1.11188e-08     -1.66e-14
polygon 363 (hole)         4 -2.60006e-08     -3.88e-14
polygon 364 (hole)        28 -1.79082e-07     -2.67e-13
polygon 365 (hole)         4 -1.58213e-10     -2.36e-16
polygon 366 (hole)         4 -9.07739e-12     -1.36e-17
polygon 367 (hole)         4 -6.73612e-09     -1.01e-14
polygon 368               16  4.14093e-01      6.18e-07
polygon 369 (hole)         4 -5.36782e-09     -8.02e-15
polygon 370 (hole)         6 -1.51450e-08     -2.26e-14
polygon 371 (hole)         4 -4.38078e-09     -6.54e-15
polygon 372 (hole)         4 -1.36409e-08     -2.04e-14
polygon 373 (hole)         4 -3.84271e-08     -5.74e-14
polygon 374 (hole)         3 -7.56649e-09     -1.13e-14
polygon 375 (hole)         4 -8.04108e-09     -1.20e-14
polygon 376 (hole)         6 -2.53481e-08     -3.78e-14
polygon 377 (hole)         4 -1.15603e-08     -1.73e-14
polygon 378 (hole)         3 -2.16577e-09     -3.23e-15
polygon 379 (hole)         3 -2.87541e-09     -4.29e-15
polygon 380 (hole)         4 -5.42479e-08     -8.10e-14
polygon 381 (hole)         3 -4.75090e-09     -7.09e-15
polygon 382 (hole)         4 -3.45047e-08     -5.15e-14
polygon 383 (hole)         6 -7.20472e-09     -1.08e-14
polygon 384 (hole)         4 -2.01300e-08     -3.01e-14
polygon 385 (hole)         9 -7.99880e-08     -1.19e-13
polygon 386 (hole)         4 -1.30938e-09     -1.96e-15
polygon 387 (hole)         9 -1.92423e-08     -2.87e-14
polygon 388 (hole)         3 -2.44412e-13     -3.65e-19
polygon 389 (hole)         4 -2.73889e-08     -4.09e-14
polygon 390 (hole)         3 -1.15109e-08     -1.72e-14
polygon 391 (hole)         4 -8.77810e-08     -1.31e-13
polygon 392               26  3.03928e+00      4.54e-06
polygon 393 (hole)         6 -8.99139e-08     -1.34e-13
polygon 394 (hole)         7 -1.69870e-07     -2.54e-13
polygon 395 (hole)        11 -1.85737e-07     -2.77e-13
polygon 396 (hole)        32 -8.16597e-07     -1.22e-12
polygon 397 (hole)         7 -6.59191e-08     -9.84e-14
polygon 398 (hole)         3 -1.48932e-08     -2.22e-14
polygon 399 (hole)        13 -3.56626e-07     -5.33e-13
polygon 400 (hole)         3 -1.64295e-12     -2.45e-18
polygon 401 (hole)         3 -2.03106e-09     -3.03e-15
polygon 402               14  1.53563e-01      2.29e-07
polygon 403 (hole)         4 -1.79041e-08     -2.67e-14
polygon 404 (hole)        15 -2.31620e-07     -3.46e-13
polygon 405 (hole)         7 -3.46670e-08     -5.18e-14
polygon 406 (hole)         6 -1.90172e-08     -2.84e-14
polygon 407 (hole)         6 -3.15828e-08     -4.72e-14
polygon 408 (hole)        13 -1.40633e-08     -2.10e-14
polygon 409 (hole)         6 -2.76240e-08     -4.12e-14
polygon 410 (hole)         3 -3.74153e-08     -5.59e-14
polygon 411 (hole)         3 -1.16808e-08     -1.74e-14
polygon 412 (hole)        10 -4.10705e-09     -6.13e-15
polygon 413 (hole)         4 -6.52777e-09     -9.75e-15
polygon 414 (hole)         6 -5.03054e-08     -7.51e-14
polygon 415 (hole)         4 -2.54303e-09     -3.80e-15
polygon 416 (hole)         3 -1.76691e-10     -2.64e-16
polygon 417 (hole)         4 -7.90017e-08     -1.18e-13
polygon 418 (hole)         3 -1.07902e-09     -1.61e-15
polygon 419 (hole)         3 -8.74115e-10     -1.31e-15
polygon 420 (hole)        12 -6.41926e-08     -9.59e-14
polygon 421 (hole)         4 -4.36895e-09     -6.52e-15
polygon 422 (hole)         4 -1.85182e-08     -2.77e-14
polygon 423 (hole)        10 -1.21248e-07     -1.81e-13
polygon 424 (hole)        10 -5.59537e-08     -8.35e-14
polygon 425 (hole)         4 -2.20079e-08     -3.29e-14
polygon 426               19  5.08538e-01      7.59e-07
polygon 427 (hole)         4 -1.45573e-08     -2.17e-14
polygon 428 (hole)         4 -4.47181e-09     -6.68e-15
polygon 429               16  1.80565e-01      2.70e-07
polygon 430 (hole)         3 -9.03958e-14     -1.35e-19
polygon 431               26  9.75091e-01      1.46e-06
polygon 432 (hole)         4 -1.90718e-09     -2.85e-15
polygon 433 (hole)         5 -1.60264e-08     -2.39e-14
polygon 434 (hole)         4 -1.82919e-08     -2.73e-14
polygon 435 (hole)         3 -7.26182e-09     -1.08e-14
polygon 436 (hole)         4 -2.12571e-08     -3.17e-14
polygon 437 (hole)         4 -2.41922e-09     -3.61e-15
polygon 438 (hole)         3 -2.80771e-08     -4.19e-14
polygon 439 (hole)         7 -5.31797e-09     -7.94e-15
polygon 440 (hole)        16 -1.29582e-07     -1.93e-13
polygon 441 (hole)         4 -1.88745e-08     -2.82e-14
polygon 442 (hole)         3 -1.45605e-08     -2.17e-14
polygon 443 (hole)         3 -2.28810e-08     -3.42e-14
polygon 444 (hole)         4 -1.59810e-08     -2.39e-14
polygon 445 (hole)         5 -1.72617e-08     -2.58e-14
polygon 446 (hole)         3 -4.96859e-09     -7.42e-15
polygon 447 (hole)        14 -8.85135e-08     -1.32e-13
polygon 448 (hole)         4 -1.85807e-08     -2.77e-14
polygon 449 (hole)         4 -8.61250e-09     -1.29e-14
polygon 450 (hole)         4 -2.06032e-08     -3.08e-14
polygon 451 (hole)         3 -1.66645e-09     -2.49e-15
polygon 452 (hole)         7 -1.66808e-08     -2.49e-14
polygon 453 (hole)         4 -1.67544e-09     -2.50e-15
polygon 454 (hole)        12 -3.97402e-08     -5.93e-14
polygon 455 (hole)         4 -4.18178e-08     -6.24e-14
polygon 456 (hole)         3 -2.93135e-09     -4.38e-15
polygon 457 (hole)         4 -6.31983e-10     -9.44e-16
polygon 458 (hole)         3 -8.50840e-13     -1.27e-18
polygon 459 (hole)         7 -8.42406e-09     -1.26e-14
polygon 460 (hole)         3 -1.34855e-12     -2.01e-18
polygon 461 (hole)         6 -1.60663e-08     -2.40e-14
polygon 462 (hole)         4 -5.66656e-11     -8.46e-17
polygon 463 (hole)         5 -1.17874e-08     -1.76e-14
polygon 464 (hole)         3 -1.75180e-14     -2.62e-20
polygon 465 (hole)         5 -2.62273e-08     -3.92e-14
polygon 466 (hole)         7 -2.12912e-08     -3.18e-14
polygon 467 (hole)        11 -4.47347e-08     -6.68e-14
polygon 468 (hole)         6 -5.00912e-09     -7.48e-15
polygon 469 (hole)         5 -3.58831e-08     -5.36e-14
polygon 470 (hole)         3 -1.06446e-09     -1.59e-15
polygon 471 (hole)         4 -1.50038e-08     -2.24e-14
polygon 472 (hole)         3 -8.74496e-10     -1.31e-15
polygon 473 (hole)         7 -3.75991e-08     -5.61e-14
polygon 474               19  2.15295e-01      3.21e-07
polygon 475 (hole)         6 -3.71337e-08     -5.54e-14
polygon 476 (hole)         3 -4.33900e-09     -6.48e-15
polygon 477 (hole)         7 -8.65067e-08     -1.29e-13
polygon 478 (hole)         8 -7.76951e-08     -1.16e-13
polygon 479 (hole)         6 -9.01907e-08     -1.35e-13
polygon 480 (hole)         8 -1.10404e-07     -1.65e-13
polygon 481 (hole)         8 -7.38263e-08     -1.10e-13
polygon 482 (hole)         3 -1.71517e-08     -2.56e-14
polygon 483               16  3.50061e-01      5.23e-07
polygon 484 (hole)         3 -1.50971e-14     -2.25e-20
polygon 485 (hole)         3 -1.78834e-04     -2.67e-10
polygon 486                3  0.00000e+00      0.00e+00
polygon 487 (hole)         4 -4.21570e-08     -6.29e-14
polygon 488 (hole)         7 -1.74178e-07     -2.60e-13
polygon 489               14  2.11571e-01      3.16e-07
polygon 490 (hole)         3 -1.05578e-08     -1.58e-14
polygon 491 (hole)         4 -4.32122e-08     -6.45e-14
polygon 492               11  1.59536e-01      2.38e-07
polygon 493 (hole)         4 -3.72418e-08     -5.56e-14
polygon 494 (hole)         6 -2.87377e-08     -4.29e-14
polygon 495 (hole)         3 -3.87472e-09     -5.79e-15
polygon 496 (hole)         6 -1.05664e-08     -1.58e-14
polygon 497 (hole)         4 -7.04519e-09     -1.05e-14
polygon 498              111  7.35028e+01      1.10e-04
polygon 499 (hole)         3 -5.99345e-09     -8.95e-15
polygon 500 (hole)         3 -2.33491e-08     -3.49e-14
polygon 501 (hole)         5 -1.14503e-07     -1.71e-13
polygon 502               19  6.34649e-01      9.48e-07
polygon 503 (hole)         3 -2.26548e-09     -3.38e-15
polygon 504 (hole)         4 -5.35016e-08     -7.99e-14
polygon 505               10  3.74827e-02      5.60e-08
polygon 506               13  1.37100e-01      2.05e-07
polygon 507               31  4.47193e+00      6.68e-06
polygon 508               19  1.23742e+00      1.85e-06
polygon 509               20  1.79201e+00      2.68e-06
polygon 510 (hole)         3 -1.64235e-08     -2.45e-14
polygon 511               16  2.19464e-01      3.28e-07
polygon 512               11  7.94092e-02      1.19e-07
polygon 513 (hole)         3 -2.99572e-08     -4.47e-14
polygon 514 (hole)         3 -8.18639e-09     -1.22e-14
polygon 515 (hole)         3 -2.81986e-08     -4.21e-14
polygon 516               14  1.45779e-01      2.18e-07
polygon 517               31  6.24672e-01      9.33e-07
polygon 518               11  2.29288e-02      3.42e-08
polygon 519               18  2.10998e-01      3.15e-07
polygon 520               26  2.52571e+00      3.77e-06
polygon 521               16  5.77783e-01      8.63e-07
polygon 522               26  8.87985e-01      1.33e-06
polygon 523               13  1.10557e-01      1.65e-07
polygon 524               21  9.55338e-01      1.43e-06
polygon 525               32  2.18189e+00      3.26e-06
polygon 526               18  7.64830e-01      1.14e-06
polygon 527               34  1.71336e+00      2.56e-06
polygon 528              115  9.29594e-01      1.39e-06
polygon 529               68  2.06761e-01      3.09e-07
polygon 530               40  4.87851e+00      7.28e-06
polygon 531               10  7.43137e-02      1.11e-07
polygon 532               17  5.95020e-01      8.88e-07
polygon 533               21  4.02132e-01      6.00e-07
polygon 534               21  1.09635e+00      1.64e-06
polygon 535               14  1.15356e-01      1.72e-07
polygon 536 (hole)        18 -1.39516e-07     -2.08e-13
polygon 537               20  3.13610e-01      4.68e-07
polygon 538              275  3.92648e+02      5.86e-04
polygon 539               16  2.19612e-01      3.28e-07
polygon 540               13  2.64858e-01      3.95e-07
polygon 541               23  1.25549e+00      1.87e-06
polygon 542                9  1.79312e-01      2.68e-07
polygon 543               15  2.93132e-01      4.38e-07
polygon 544               15  5.29892e-01      7.91e-07
polygon 545               23  4.95351e-01      7.40e-07
polygon 546 (hole)        10 -2.23148e-07     -3.33e-13
polygon 547               14  4.41976e-01      6.60e-07
polygon 548               19  3.32318e-01      4.96e-07
polygon 549 (hole)         6 -6.70431e-08     -1.00e-13
polygon 550               16  4.15914e-01      6.21e-07
polygon 551               22  5.10526e-01      7.62e-07
polygon 552               48  1.87148e+01      2.79e-05
polygon 553               18  1.62882e+00      2.43e-06
polygon 554 (hole)         9 -1.96623e-07     -2.94e-13
polygon 555               71  2.43615e+01      3.64e-05
polygon 556              116  5.16444e+01      7.71e-05
polygon 557               19  1.54144e+00      2.30e-06
polygon 558               11  2.01013e-01      3.00e-07
polygon 559 (hole)         4 -1.87899e-08     -2.81e-14
polygon 560               20  2.31201e+00      3.45e-06
polygon 561               15  5.94674e-01      8.88e-07
polygon 562              218  1.56252e+02      2.33e-04
polygon 563              162  8.80247e+01      1.31e-04
polygon 564               16  3.69968e-01      5.52e-07
polygon 565               14  1.61451e-01      2.41e-07
polygon 566               18  5.54566e-01      8.28e-07
polygon 567               19  5.07078e-01      7.57e-07
polygon 568 (hole)         7 -1.85276e-07     -2.77e-13
polygon 569               57  6.86475e+00      1.03e-05
polygon 570               15  6.56018e-01      9.80e-07
polygon 571               21  6.65092e-01      9.93e-07
polygon 572               15  5.97784e-01      8.93e-07
polygon 573               37  2.56010e+00      3.82e-06
polygon 574               13  1.66996e-01      2.49e-07
polygon 575               20  8.18381e-01      1.22e-06
polygon 576               21  2.89704e+00      4.33e-06
polygon 577               15  4.08213e-01      6.10e-07
polygon 578 (hole)         3 -1.98618e-09     -2.97e-15
polygon 579              104  4.71547e+01      7.04e-05
polygon 580               28  3.80443e+00      5.68e-06
polygon 581               67  3.17343e+01      4.74e-05
polygon 582              730  9.20735e+02      1.37e-03
polygon 583               16  2.01576e-01      3.01e-07
polygon 584               12  1.41506e-01      2.11e-07
polygon 585               15  4.92125e-01      7.35e-07
polygon 586               18  1.58128e+00      2.36e-06
polygon 587 (hole)         7 -6.67005e-08     -9.96e-14
polygon 588               25  1.94322e+00      2.90e-06
polygon 589               12  2.51373e-01      3.75e-07
polygon 590              412  4.47936e+02      6.69e-04
polygon 591               22  7.57105e-01      1.13e-06
polygon 592               16  6.97434e-01      1.04e-06
polygon 593 (hole)         3 -2.86672e-08     -4.28e-14
polygon 594               39  4.45130e+00      6.65e-06
polygon 595 (hole)         4 -2.17002e-08     -3.24e-14
polygon 596 (hole)         3 -1.55741e-08     -2.33e-14
polygon 597               13  8.83786e-01      1.32e-06
polygon 598               10  2.06200e-01      3.08e-07
polygon 599 (hole)         3 -3.50430e-08     -5.23e-14
polygon 600               76  3.60497e+01      5.38e-05
polygon 601 (hole)         5 -5.41615e-08     -8.09e-14
polygon 602               13  4.56433e-01      6.82e-07
polygon 603               28  2.30613e+00      3.44e-06
polygon 604 (hole)         4 -1.59313e-08     -2.38e-14
polygon 605               27  2.50338e+00      3.74e-06
polygon 606               14  4.69238e-01      7.01e-07
polygon 607               10  1.69886e-01      2.54e-07
polygon 608 (hole)         8 -1.09145e-07     -1.63e-13
polygon 609               16  5.42822e-01      8.11e-07
polygon 610               10  1.17633e-01      1.76e-07
polygon 611                8  7.08579e-02      1.06e-07
polygon 612              144  7.85300e+01      1.17e-04
polygon 613               37  1.93477e+00      2.89e-06
polygon 614 (hole)         3 -1.47049e-09     -2.20e-15
polygon 615               25  5.59996e-01      8.36e-07
polygon 616 (hole)         8 -1.20535e-07     -1.80e-13
polygon 617               44  2.86031e+00      4.27e-06
polygon 618               26  1.26276e+00      1.89e-06
polygon 619              149  1.37840e+02      2.06e-04
polygon 620               18  4.84958e-01      7.24e-07
polygon 621               73  2.93195e+01      4.38e-05
polygon 622 (hole)         5 -1.99458e-08     -2.98e-14
polygon 623               45  6.87481e+00      1.03e-05
polygon 624 (hole)         4 -1.81549e-08     -2.71e-14
polygon 625               64  1.80880e+01      2.70e-05
polygon 626               12  8.76879e-01      1.31e-06
polygon 627               26  2.78381e+00      4.16e-06
polygon 628               35  8.72326e+00      1.30e-05
polygon 629               53  1.06237e+01      1.59e-05
polygon 630               26  5.40467e+00      8.07e-06
polygon 631              148  1.05037e+02      1.57e-04
polygon 632               21  8.32798e-01      1.24e-06
polygon 633               21  2.23023e+00      3.33e-06
polygon 634                8  6.30805e-01      9.42e-07
polygon 635               78  3.67603e+01      5.49e-05
polygon 636               18  1.64745e+00      2.46e-06
polygon 637               58  1.31747e+01      1.97e-05
polygon 638               94  1.16837e+01      1.74e-05
polygon 639               53  3.17801e+00      4.75e-06
polygon 640              137  8.46454e+00      1.26e-05
polygon 641               23  4.36507e-01      6.52e-07
polygon 642               27  9.72136e-01      1.45e-06
polygon 643               55  1.98854e+00      2.97e-06
polygon 644               48  1.02651e+01      1.53e-05
polygon 645               15  3.16153e-01      4.72e-07
polygon 646               22  8.93201e-01      1.33e-06
polygon 647               18  2.91538e-01      4.35e-07
polygon 648               59  3.70736e+00      5.54e-06
polygon 649               11  1.55205e-01      2.32e-07
polygon 650               14  1.87401e-01      2.80e-07
polygon 651               11  6.54128e-02      9.77e-08
polygon 652               11  8.64659e-02      1.29e-07
polygon 653               12  2.25129e-01      3.36e-07
polygon 654               30  2.72178e+00      4.06e-06
polygon 655 (hole)         3 -2.26252e-06     -3.38e-12
polygon 656              180  1.80314e+01      2.69e-05
polygon 657               14  4.01942e-01      6.00e-07
polygon 658              103  1.33467e+01      1.99e-05
polygon 659               14  2.07708e-01      3.10e-07
polygon 660               55  4.38623e+00      6.55e-06
polygon 661               41  4.69733e+00      7.01e-06
polygon 662               74  3.92633e+00      5.86e-06
polygon 663               31  2.79477e+00      4.17e-06
polygon 664               17  8.74607e-01      1.31e-06
polygon 665               26  4.50962e-01      6.73e-07
polygon 666               18  1.15065e+00      1.72e-06
polygon 667              164  1.48037e+01      2.21e-05
polygon 668               12  3.17290e-01      4.74e-07
polygon 669               69  2.59551e+01      3.88e-05
polygon 670               37  7.85963e-01      1.17e-06
polygon 671               14  3.37873e-01      5.05e-07
polygon 672               12  1.77980e-01      2.66e-07
polygon 673               21  6.46195e-01      9.65e-07
polygon 674               24  2.69998e+00      4.03e-06
polygon 675               30  1.00675e+00      1.50e-06
polygon 676               57  4.22747e+00      6.31e-06
polygon 677               23  8.63543e-01      1.29e-06
polygon 678               11  2.01925e-01      3.02e-07
polygon 679               95  2.88232e+01      4.30e-05
polygon 680              452  9.64208e+01      1.44e-04
polygon 681               77  1.72297e+01      2.57e-05
polygon 682               16  3.33522e-01      4.98e-07
polygon 683               19  6.82895e-01      1.02e-06
polygon 684               43  3.95062e+00      5.90e-06
polygon 685               29  1.00709e+00      1.50e-06
polygon 686               21  8.90014e-01      1.33e-06
polygon 687               22  7.93453e-01      1.18e-06
polygon 688               13  4.52813e-01      6.76e-07
polygon 689               30  1.49817e+00      2.24e-06
polygon 690               48  4.23996e+00      6.33e-06
polygon 691               37  1.20290e+00      1.80e-06
polygon 692               19  4.85805e-01      7.25e-07
polygon 693               46  2.18001e+00      3.26e-06
polygon 694                7  1.33723e-01      2.00e-07
polygon 695               54  1.82059e+01      2.72e-05
polygon 696               10  8.57866e-01      1.28e-06
polygon 697               14  1.98445e-01      2.96e-07
polygon 698               19  4.68357e-01      6.99e-07
polygon 699               52  8.57765e+00      1.28e-05
polygon 700               73  7.91998e+00      1.18e-05
polygon 701               23  8.19561e+00      1.22e-05
polygon 702               10  1.06594e-01      1.59e-07
polygon 703              169  1.39462e+01      2.08e-05
polygon 704               17  2.55915e-01      3.82e-07
polygon 705               12  1.81516e-01      2.71e-07
polygon 706               16  3.68509e-01      5.50e-07
polygon 707               12  1.24954e-01      1.87e-07
polygon 708             1012  4.61886e+02      6.90e-04
polygon 709               15  5.42647e-01      8.10e-07
polygon 710               12  5.53088e-01      8.26e-07
polygon 711               24  1.06455e+00      1.59e-06
polygon 712               14  3.17038e-01      4.73e-07
polygon 713               17  1.25559e+00      1.87e-06
polygon 714               16  1.92407e+00      2.87e-06
polygon 715               80  1.86306e+01      2.78e-05
polygon 716               11  6.83880e-02      1.02e-07
polygon 717              166  4.21470e+01      6.29e-05
polygon 718               22  1.36174e+00      2.03e-06
polygon 719               31  6.07363e-01      9.07e-07
polygon 720               10  1.27067e-01      1.90e-07
polygon 721              164  5.02951e+01      7.51e-05
polygon 722               17  1.20798e+00      1.80e-06
polygon 723               59  2.41720e+00      3.61e-06
polygon 724               12  2.52952e-01      3.78e-07
polygon 725               14  2.44407e-01      3.65e-07
polygon 726               43  4.89177e+00      7.30e-06
polygon 727               40  3.86642e+00      5.77e-06
polygon 728               23  4.59203e-01      6.86e-07
polygon 729               14  1.79442e-01      2.68e-07
polygon 730               13  1.61191e-01      2.41e-07
polygon 731               43  1.33365e+00      1.99e-06
polygon 732               75  1.20247e+01      1.80e-05
polygon 733               17  9.96466e-02      1.49e-07
polygon 734               22  1.60660e+00      2.40e-06
polygon 735              730  1.16983e+02      1.75e-04
polygon 736              129  7.87539e+00      1.18e-05
polygon 737               22  4.66256e-01      6.96e-07
polygon 738               19  5.44519e-01      8.13e-07
polygon 739               16  3.73911e-01      5.58e-07
polygon 740               88  3.47678e+01      5.19e-05
polygon 741               46  3.23711e+00      4.83e-06
polygon 742               47  9.17482e-01      1.37e-06
polygon 743               14  4.44869e-01      6.64e-07
polygon 744               43  1.29527e+00      1.93e-06
polygon 745               65  3.20974e+00      4.79e-06
polygon 746              306  1.83695e+01      2.74e-05
polygon 747               17  4.25486e-01      6.35e-07
polygon 748               29  1.08467e+00      1.62e-06
polygon 749               26  1.16602e+00      1.74e-06
polygon 750              135  4.90655e+00      7.33e-06
polygon 751               23  1.87723e+00      2.80e-06
polygon 752               21  7.02338e-01      1.05e-06
polygon 753               30  4.99970e+00      7.47e-06
polygon 754               21  7.71571e-01      1.15e-06
polygon 755               13  1.47832e-01      2.21e-07
polygon 756               26  3.17307e+00      4.74e-06
polygon 757               17  1.03642e+00      1.55e-06
polygon 758              100  8.74537e+00      1.31e-05
polygon 759               11  1.13950e-01      1.70e-07
polygon 760               26  8.27840e-01      1.24e-06
polygon 761               17  4.10882e-01      6.14e-07
polygon 762               13  1.49250e-01      2.23e-07
polygon 763               35  1.76408e+00      2.63e-06
polygon 764               15  2.03390e-01      3.04e-07
polygon 765               18  1.70982e-01      2.55e-07
polygon 766               32  7.98599e+00      1.19e-05
polygon 767               27  3.20839e-01      4.79e-07
polygon 768              539  1.17139e+02      1.75e-04
polygon 769               51  1.57614e+00      2.35e-06
polygon 770               19  3.60772e-01      5.39e-07
polygon 771               58  2.10060e+00      3.14e-06
polygon 772               13  1.70224e-01      2.54e-07
polygon 773               11  6.11736e-01      9.13e-07
polygon 774               35  5.26152e+00      7.86e-06
polygon 775               16  1.99371e-01      2.98e-07
polygon 776               14  4.48979e-01      6.70e-07
polygon 777               25  1.92913e+00      2.88e-06
polygon 778               16  4.26157e-01      6.36e-07
polygon 779               46  1.88229e+00      2.81e-06
polygon 780               17  2.29456e-01      3.43e-07
polygon 781               36  4.01418e+00      5.99e-06
polygon 782               57  1.77153e+01      2.65e-05
polygon 783               20  5.72298e-01      8.55e-07
polygon 784               14  2.46782e-01      3.68e-07
polygon 785               29  1.01300e+00      1.51e-06
polygon 786               90  8.54955e+00      1.28e-05
polygon 787               21  3.99098e-01      5.96e-07
polygon 788               12  2.81650e-01      4.21e-07
polygon 789               52  3.19337e+00      4.77e-06
polygon 790               13  3.40300e-01      5.08e-07
polygon 791               29  9.08533e-01      1.36e-06
polygon 792               20  3.07393e-01      4.59e-07
polygon 793               19  1.15906e+00      1.73e-06
polygon 794               17  1.16182e+00      1.73e-06
polygon 795               21  2.60877e+00      3.90e-06
polygon 796               15  1.09001e-01      1.63e-07
polygon 797               37  1.31921e+00      1.97e-06
polygon 798               53  1.46854e+00      2.19e-06
polygon 799               59  8.92822e+00      1.33e-05
polygon 800                9  7.44981e-02      1.11e-07
polygon 801               14  2.41462e-01      3.61e-07
polygon 802               96  7.25940e+00      1.08e-05
polygon 803               11  1.06055e-01      1.58e-07
polygon 804               49  1.87834e+00      2.80e-06
polygon 805               23  6.08310e-01      9.08e-07
polygon 806               50  6.76488e-01      1.01e-06
polygon 807               22  4.83089e-01      7.21e-07
polygon 808               17  1.17278e-01      1.75e-07
polygon 809               13  8.98786e-01      1.34e-06
polygon 810               43  1.01757e+00      1.52e-06
polygon 811               52  1.68377e+00      2.51e-06
polygon 812              348  2.50314e+02      3.74e-04
polygon 813               43  1.29120e+00      1.93e-06
polygon 814               71  2.02836e+00      3.03e-06
polygon 815 (hole)         4 -4.12791e-06     -6.16e-12
polygon 816              141  4.53240e+00      6.77e-06
polygon 817               52  3.62008e+00      5.41e-06
polygon 818               20  7.69539e-01      1.15e-06
polygon 819               88  8.88904e+00      1.33e-05
polygon 820               12  1.44668e-01      2.16e-07
polygon 821               46  2.28026e+00      3.40e-06
polygon 822               39  4.79165e+00      7.15e-06
polygon 823               53  2.86736e+00      4.28e-06
polygon 824               38  3.22508e+00      4.82e-06
polygon 825               18  4.44863e-01      6.64e-07
polygon 826               45  1.45134e+00      2.17e-06
polygon 827               86  2.56400e+01      3.83e-05
polygon 828              158  1.74510e+01      2.61e-05
polygon 829               14  2.50383e-01      3.74e-07
polygon 830              120  7.17019e+00      1.07e-05
polygon 831              146  2.28448e+01      3.41e-05
polygon 832              143  2.34472e+01      3.50e-05
polygon 833               20  3.47415e-01      5.19e-07
polygon 834 (hole)         4 -4.70080e-12     -7.02e-18
polygon 835               84  3.86387e+00      5.77e-06
polygon 836               68  1.10850e+01      1.66e-05
polygon 837              863  7.47702e+01      1.12e-04
polygon 838               53  7.54742e+00      1.13e-05
polygon 839 (hole)         5 -5.07266e-12     -7.57e-18
polygon 840               64  2.11016e+00      3.15e-06
polygon 841              103  1.47525e+01      2.20e-05
polygon 842               26  3.40345e+00      5.08e-06
polygon 843 (hole)         6 -3.63322e-12     -5.43e-18
polygon 844 (hole)         8 -7.62695e-12     -1.14e-17
polygon 845              103  8.40445e+00      1.25e-05
polygon 846               23  4.92041e-01      7.35e-07
polygon 847               27  5.53625e-01      8.27e-07
polygon 848               60  3.56102e+00      5.32e-06
polygon 849              151  1.51314e+01      2.26e-05
polygon 850 (hole)         3 -7.39478e-13     -1.10e-18
polygon 851              708  4.18129e+02      6.24e-04
polygon 852              120  1.32857e+01      1.98e-05
polygon 853              856  2.53898e+02      3.79e-04
polygon 854 (hole)         3 -1.51167e-12     -2.26e-18
polygon 855 (hole)         3 -4.16327e-13     -6.22e-19
polygon 856 (hole)         3 -1.81968e-13     -2.72e-19
polygon 857 (hole)         5 -8.69542e-12     -1.30e-17
polygon 858               13  2.47298e-01      3.69e-07
polygon 859 (hole)         5 -7.23071e-12     -1.08e-17
polygon 860               16  7.35754e-01      1.10e-06
polygon 861 (hole)         4 -5.36614e-12     -8.01e-18
polygon 862               21  5.31633e-01      7.94e-07
polygon 863 (hole)        11 -4.03802e-12     -6.03e-18
polygon 864 (hole)         8 -8.21058e-12     -1.23e-17
polygon 865 (hole)         3 -8.94235e-13     -1.34e-18
polygon 866 (hole)         5 -8.77891e-13     -1.31e-18
polygon 867               72  1.03185e+01      1.54e-05
polygon 868 (hole)         4 -1.46473e-12     -2.19e-18
polygon 869 (hole)         3 -9.57738e-14     -1.43e-19
polygon 870               38  1.34734e+00      2.01e-06
polygon 871                8  4.09489e-01      6.11e-07
polygon 872 (hole)         4 -2.72557e-12     -4.07e-18
polygon 873               17  2.46995e-01      3.69e-07
polygon 874 (hole)         3 -4.56587e-13     -6.82e-19
polygon 875               14  1.42285e-01      2.12e-07
polygon 876 (hole)         4 -1.74558e-12     -2.61e-18
polygon 877 (hole)         6 -5.91253e-12     -8.83e-18
polygon 878 (hole)         4 -1.34510e-13     -2.01e-19
polygon 879 (hole)         3 -6.12765e-17     -9.15e-23
polygon 880                9  1.56346e-01      2.33e-07
polygon 881                4  2.90958e-17      4.34e-23
polygon 882               42  5.49581e-01      8.21e-07
polygon 883 (hole)         5 -2.12744e-12     -3.18e-18
polygon 884 (hole)         3 -5.23582e-12     -7.82e-18
polygon 885 (hole)         3 -2.43443e-13     -3.64e-19
polygon 886 (hole)         3 -8.02710e-13     -1.20e-18
polygon 887 (hole)         6 -1.09890e-11     -1.64e-17
polygon 888 (hole)         4 -1.53858e-13     -2.30e-19
polygon 889               59  1.15813e+00      1.73e-06
polygon 890 (hole)         4 -2.30355e-12     -3.44e-18
polygon 891               12  1.16845e-01      1.74e-07
polygon 892               11  6.00122e-02      8.96e-08
polygon 893               14  1.63465e-01      2.44e-07
polygon 894                8  4.73674e-02      7.07e-08
polygon 895               23  6.69330e-01      9.99e-07
polygon 896 (hole)         3 -1.13772e-13     -1.70e-19
polygon 897              144  6.48346e+00      9.68e-06
polygon 898               74  4.79923e+00      7.17e-06
polygon 899              116  1.02706e+01      1.53e-05
polygon 900               16  2.74005e-01      4.09e-07
polygon 901               97  4.72599e+00      7.06e-06
polygon 902               18  6.15583e-01      9.19e-07
polygon 903              146  1.38979e+01      2.08e-05
polygon 904 (hole)         5 -4.08205e-12     -6.10e-18
polygon 905               32  1.94516e+00      2.90e-06
polygon 906               26  9.81356e-01      1.47e-06
polygon 907 (hole)         3 -2.01834e-12     -3.01e-18
polygon 908 (hole)         4 -9.94985e-13     -1.49e-18
polygon 909               82  4.39797e+00      6.57e-06
polygon 910               33  9.36113e-01      1.40e-06
polygon 911 (hole)         3 -1.40697e-13     -2.10e-19
polygon 912               47  1.51844e+00      2.27e-06
polygon 913 (hole)         3 -6.88122e-12     -1.03e-17
polygon 914              196  2.32689e+01      3.47e-05
polygon 915 (hole)         3 -2.40592e-12     -3.59e-18
polygon 916 (hole)         5 -6.19998e-12     -9.26e-18
polygon 917               53  1.30069e+00      1.94e-06
polygon 918               70  1.00802e+01      1.51e-05
polygon 919 (hole)         3 -3.17405e-13     -4.74e-19
polygon 920               10  5.33870e-01      7.97e-07
polygon 921              211  2.73026e+01      4.08e-05
polygon 922               41  1.44141e+00      2.15e-06
polygon 923                9  5.20220e-01      7.77e-07
polygon 924              316  6.66808e+01      9.96e-05
polygon 925               17  7.48505e-01      1.12e-06
polygon 926               15  1.49799e-01      2.24e-07
polygon 927               24  4.99500e-01      7.46e-07
polygon 928               12  2.24463e-01      3.35e-07
polygon 929              100  1.06028e+01      1.58e-05
polygon 930               11  1.24267e-01      1.86e-07
polygon 931              143  1.89581e+01      2.83e-05
polygon 932               60  7.07802e+00      1.06e-05
polygon 933              124  1.23476e+01      1.84e-05
polygon 934               34  5.67625e+00      8.48e-06
polygon 935              797  1.86955e+02      2.79e-04
polygon 936               15  7.15891e-01      1.07e-06
polygon 937               20  1.61624e+00      2.41e-06
polygon 938               13  1.73955e-01      2.60e-07
polygon 939               60  2.42590e+00      3.62e-06
polygon 940               20  3.86365e-01      5.77e-07
polygon 941              213  3.69421e+01      5.52e-05
polygon 942               15  4.55574e-01      6.80e-07
polygon 943               35  1.56058e+00      2.33e-06
polygon 944               32  1.42419e+00      2.13e-06
polygon 945               68  2.30711e+00      3.44e-06
polygon 946              240  5.93935e+01      8.87e-05
polygon 947              145  1.58875e+01      2.37e-05
polygon 948               39  3.67239e+00      5.48e-06
polygon 949               54  3.84847e+00      5.75e-06
polygon 950               32  1.14073e+00      1.70e-06
polygon 951               18  3.58549e-01      5.35e-07
polygon 952               25  8.06470e-01      1.20e-06
polygon 953               15  2.67818e-01      4.00e-07
polygon 954               26  1.43681e+00      2.15e-06
polygon 955               18  3.10061e-01      4.63e-07
polygon 956               18  4.71644e-01      7.04e-07
polygon 957               16  5.04146e-01      7.53e-07
polygon 958               72  7.33720e+00      1.10e-05
polygon 959               27  1.37772e+00      2.06e-06
polygon 960               15  3.17217e-01      4.74e-07
polygon 961               37  1.17498e+00      1.75e-06
polygon 962               21  7.05388e-01      1.05e-06
polygon 963              216  2.05399e+01      3.07e-05
polygon 964               16  2.40093e-01      3.59e-07
polygon 965               29  1.71282e+00      2.56e-06
polygon 966               30  1.28695e+00      1.92e-06
polygon 967               36  2.29670e+00      3.43e-06
polygon 968              129  1.97698e+01      2.95e-05
polygon 969               65  2.65969e+00      3.97e-06
polygon 970              283  3.30575e+01      4.94e-05
polygon 971               28  7.06271e-01      1.05e-06
polygon 972               28  5.50158e-01      8.21e-07
polygon 973               49  2.25015e+00      3.36e-06
polygon 974               26  1.24280e+00      1.86e-06
polygon 975               28  9.02794e-01      1.35e-06
polygon 976               18  4.84064e-01      7.23e-07
polygon 977               47  2.33442e+00      3.49e-06
polygon 978               17  2.23835e-01      3.34e-07
polygon 979               34  7.32504e-01      1.09e-06
polygon 980               41  7.15538e-01      1.07e-06
polygon 981               21  5.81687e-01      8.69e-07
polygon 982               63  3.83819e+00      5.73e-06
polygon 983               27  1.33192e+00      1.99e-06
polygon 984               46  2.95028e+00      4.41e-06
polygon 985               10  1.34210e-01      2.00e-07
polygon 986               16  2.38274e-01      3.56e-07
polygon 987               18  4.18536e-01      6.25e-07
polygon 988               46  1.26584e+00      1.89e-06
polygon 989               14  2.14679e-01      3.21e-07
polygon 990               76  3.63371e+00      5.43e-06
polygon 991              339  4.44685e+01      6.64e-05
polygon 992              109  6.29786e+00      9.40e-06
polygon 993               14  2.89570e-01      4.32e-07
polygon 994               37  9.07704e-01      1.36e-06
polygon 995               68  3.86104e+00      5.77e-06
polygon 996              183  1.45458e+01      2.17e-05
polygon 997               44  1.64355e+00      2.45e-06
polygon 998               26  1.21667e+00      1.82e-06
polygon 999               13  1.95710e-01      2.92e-07
polygon 1000              38  1.44402e+00      2.16e-06
polygon 1001              14  3.29394e-01      4.92e-07
polygon 1002              14  1.79828e-01      2.69e-07
polygon 1003              44  2.11118e+00      3.15e-06
polygon 1004             643  1.79215e+02      2.68e-04
polygon 1005              24  1.10979e+00      1.66e-06
polygon 1006              30  1.57106e+00      2.35e-06
polygon 1007              46  3.61716e+00      5.40e-06
polygon 1008              93  7.56340e-01      1.13e-06
polygon 1009              57  3.07732e+00      4.59e-06
polygon 1010              18  2.16224e-01      3.23e-07
polygon 1011              17  5.12388e-01      7.65e-07
polygon 1012               9  7.44780e-02      1.11e-07
polygon 1013              57  8.09122e+00      1.21e-05
polygon 1014              16  4.01323e-01      5.99e-07
polygon 1015             118  2.90304e+01      4.33e-05
polygon 1016              22  4.44693e-01      6.64e-07
polygon 1017              43  2.28763e+00      3.42e-06
polygon 1018              27  1.15444e+00      1.72e-06
polygon 1019              25  1.21695e+00      1.82e-06
polygon 1020              53  4.49228e+00      6.71e-06
polygon 1021 (hole)        8 -5.14959e-08     -7.69e-14
polygon 1022 (hole)        3 -1.87101e-12     -2.79e-18
polygon 1023 (hole)       31 -2.12239e-07     -3.17e-13
polygon 1024 (hole)        6 -4.99400e-08     -7.46e-14
polygon 1025 (hole)        8 -1.14003e-07     -1.70e-13
polygon 1026 (hole)        4 -5.05161e-08     -7.54e-14
polygon 1027 (hole)        4 -6.87492e-09     -1.03e-14
polygon 1028 (hole)        6 -5.81106e-12     -8.68e-18
polygon 1029 (hole)        7 -1.10596e-07     -1.65e-13
polygon 1030 (hole)       14 -3.52603e-07     -5.26e-13
polygon 1031 (hole)       14 -1.57740e-07     -2.36e-13
polygon 1032 (hole)        6 -6.42551e-08     -9.59e-14
polygon 1033 (hole)        4 -3.16857e-08     -4.73e-14
polygon 1034 (hole)        8 -7.25950e-08     -1.08e-13
polygon 1035 (hole)        4 -1.10823e-08     -1.65e-14
polygon 1036 (hole)        8 -1.54074e-07     -2.30e-13
polygon 1037 (hole)        4 -5.15531e-09     -7.70e-15
polygon 1038 (hole)        7 -4.08919e-08     -6.11e-14
polygon 1039 (hole)        9 -6.45354e-08     -9.64e-14
polygon 1040 (hole)       13 -6.64876e-08     -9.93e-14
polygon 1041 (hole)        5 -5.31928e-08     -7.94e-14
polygon 1042 (hole)        6 -4.67379e-08     -6.98e-14
polygon 1043 (hole)        4 -2.58840e-08     -3.86e-14
polygon 1044 (hole)       10 -9.35297e-08     -1.40e-13
polygon 1045 (hole)        4 -1.76462e-08     -2.63e-14
polygon 1046 (hole)        4 -5.14884e-08     -7.69e-14
polygon 1047 (hole)        3 -1.43941e-09     -2.15e-15
polygon 1048 (hole)        8 -2.11864e-08     -3.16e-14
polygon 1049 (hole)        4 -3.97950e-09     -5.94e-15
polygon 1050 (hole)        3 -2.72931e-08     -4.08e-14
polygon 1051 (hole)        4 -9.91931e-08     -1.48e-13
polygon 1052 (hole)        9 -4.30090e-08     -6.42e-14
polygon 1053 (hole)        3 -2.69442e-08     -4.02e-14
polygon 1054 (hole)        6 -3.65336e-08     -5.46e-14
polygon 1055 (hole)       10 -5.24779e-08     -7.84e-14
polygon 1056 (hole)       12 -5.88191e-08     -8.78e-14
polygon 1057 (hole)        8 -3.34285e-08     -4.99e-14
polygon 1058 (hole)       18 -1.75177e-07     -2.62e-13
polygon 1059 (hole)        6 -6.63687e-08     -9.91e-14
polygon 1060 (hole)        3 -1.24000e-08     -1.85e-14
polygon 1061 (hole)       10 -1.14560e-07     -1.71e-13
polygon 1062 (hole)        5 -5.13944e-08     -7.67e-14
polygon 1063 (hole)       10 -8.30360e-08     -1.24e-13
polygon 1064 (hole)        4 -1.89513e-08     -2.83e-14
polygon 1065 (hole)        6 -1.07333e-07     -1.60e-13
polygon 1066 (hole)        3 -5.80811e-09     -8.67e-15
polygon 1067              79  1.47390e+01      2.20e-05
polygon 1068              54  9.35305e+00      1.40e-05
polygon 1069 (hole)        3 -1.94488e-09     -2.90e-15
polygon 1070 (hole)       10 -6.90517e-08     -1.03e-13
polygon 1071 (hole)        4 -1.17368e-08     -1.75e-14
polygon 1072 (hole)        9 -5.79029e-08     -8.65e-14
polygon 1073 (hole)        3 -4.78872e-08     -7.15e-14
polygon 1074 (hole)        3 -5.90058e-09     -8.81e-15
polygon 1075 (hole)        5 -6.11108e-08     -9.12e-14
polygon 1076           37606  6.60254e+05      9.86e-01
polygon 1077 (hole)        5 -2.16281e-12     -3.23e-18
enclosing rectangle: [-210.0086, 724.6476] x [1072.0263, 3158.4671] km
                     (934.7 x 2086 km)
Window area = 669714 square km
Unit of length: 1 km
Fraction of frame area: 0.343

The ppp object outputted from combining both the point and polygon feature results in the boundary of Myanmar outlining the plot of conflict events as shown.

# Plot each masked ppp object
par(mfrow = c(2,3), mar = c(0,0,1,0))
for (quarter in names(masked_ppp_list_km)) {
  plot(masked_ppp_list_km[[quarter]], main = paste("Year Quarter:", quarter))
}

Let’s also create combined ppp and owin objects by the districts I am interested in this study. We’ll prepare the objects for the districts: Yinmarbin, Shwebo, Pakokku and Mandalay.

# Prepare Dataset for Yinmarbin District
conflict_yinmarbin = filter(conflict_data_sf, DT == "Yinmarbin")
boundary_yinmarbin <- filter(boundary_sf, DT == "Yinmarbin")

# Create a combined ppp and owin object
yinmarbin_owin <- as.owin(boundary_yinmarbin)
ppp_obj <- as.ppp(conflict_yinmarbin$geometry)

ppp_obj <- as.ppp(st_geometry(conflict_yinmarbin))

# Handle duplicates
ppp_obj <- rjitter(ppp_obj, retry=TRUE, nsim=1, drop=TRUE)

yinmarbin_ppp_owin <- ppp_obj[yinmarbin_owin]
yinmarbin_ppp_owin <- rescale(yinmarbin_ppp_owin, 1000, "km")
# Prepare Dataset for Shwebo District
conflict_shwebo = filter(conflict_data_sf, DT == "Shwebo")
boundary_shwebo <- filter(boundary_sf, DT == "Shwebo")

# Create a combined ppp and owin object
shwebo_owin <- as.owin(boundary_shwebo)
ppp_obj <- as.ppp(conflict_shwebo$geometry)

# Handle duplicates
ppp_obj <- rjitter(ppp_obj, retry=TRUE, nsim=1, drop=TRUE)

shwebo_ppp_owin <- ppp_obj[shwebo_owin]
shwebo_ppp_owin <- rescale(shwebo_ppp_owin, 1000, "km")
# Prepare Dataset for Pakokku District
conflict_pakokku = filter(conflict_data_sf, DT == "Pakokku")
boundary_pakokku <- filter(boundary_sf, DT == "Pakokku")

# Create a combined ppp and owin object
pakokku_owin <- as.owin(boundary_pakokku)
ppp_obj <- as.ppp(conflict_pakokku$geometry)

# Handle duplicates
ppp_obj <- rjitter(ppp_obj, retry=TRUE, nsim=1, drop=TRUE)

pakokku_ppp_owin <- ppp_obj[pakokku_owin]
pakokku_ppp_owin <- rescale(pakokku_ppp_owin, 1000, "km")
# Prepare Dataset for Mandalay District
conflict_mandalay = filter(conflict_data_sf, DT == "Mandalay")
boundary_mandalay <- filter(boundary_sf, DT == "Mandalay")

# Create a combined ppp and owin object
mandalay_owin <- as.owin(boundary_mandalay)
ppp_obj <- as.ppp(conflict_mandalay$geometry)

# Handle duplicates
ppp_obj <- rjitter(ppp_obj, retry=TRUE, nsim=1, drop=TRUE)

mandalay_ppp_owin <- ppp_obj[mandalay_owin]
mandalay_ppp_owin <- rescale(mandalay_ppp_owin, 1000, "km")
par(mfrow = c(2,2), mar = c(0,0,1,0))
plot(yinmarbin_ppp_owin, main = paste("Yinmarbin District"))
plot(shwebo_ppp_owin, main = paste("Shwebo District"))
plot(pakokku_ppp_owin, main = paste("Pakokku District"))
plot(mandalay_ppp_owin, main = paste("Mandalay District"))

Observations

Here I plot the districts with the highest concentration of armed conflicts. We can observe localised clustering of armed conflicts in certain towns of each district (darker spots on the map) which also indicates high frequency of conflicts occurring in these areas. On the other hand, there are also parts of each district with smaller frequency of conflicts but clusters of conflicts formed.

4. Exploratory Data Analysis

4.1 Identifying Districts with Highest Proportion of Conflicts

It’ll also be interesting to find out specific districts with the highest concentration of armed conflicts. I will first calculate the total occurrences of conflict events per district and add the column to boundary_sf.

Count number of conflicts by districts
conflict_count <- conflict_data_sf %>%
  group_by(DT) %>%
  summarise(total_count_DT = n()) %>%
  st_drop_geometry() %>%
  select(DT, total_count_DT)

# Perform the join
boundary_sf <- boundary_sf %>%
  left_join(conflict_count, by = "DT")

Next, let’s calculate the proportion of total conflicts and add it as a column into the boundary_sf dataset as proportion_DT.

# Create new 'proportion_DT' column
boundary_sf <- boundary_sf %>%
  mutate(proportion_DT = total_count_DT / sum(total_count_DT))
head(boundary_sf[c('DT','total_count_DT','proportion_DT')])
Simple feature collection with 6 features and 3 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -14915.04 ymin: 1736124 xmax: 187961.7 ymax: 2051144
Projected CRS: WGS 84 / UTM zone 47N
         DT total_count_DT proportion_DT                       geometry
1  Hinthada            160            NA MULTIPOLYGON (((90859.89 20...
2   Labutta             51            NA MULTIPOLYGON (((75991.51 17...
3    Maubin            118            NA MULTIPOLYGON (((115559 1928...
4 Myaungmya             59            NA MULTIPOLYGON (((39919.39 18...
5   Pathein            334            NA MULTIPOLYGON (((-6302.348 1...
6    Pyapon            131            NA MULTIPOLYGON (((93411.72 17...

At a quick glance, we can see that central and southern parts of Myanmar have the highest proportions of armed conflict events occurring.

Set up the points map
districts_choropleth <-
tm_shape(boundary_sf) +
  tm_fill("proportion_DT",
          n=10,
          title="Proportion",
          style="equal",
          palette="Blues") +
  tm_borders(lwd=0.2,
             alpha=1) +
  tm_text(text = "DT", 
          size = 0.2, 
          col = "black",
          fontface = "bold") +
  tm_layout(main.title = "Distribution of Conflict Points Across Districts",
            legend.outside=FALSE,
            main.title.size=1)
# Plot the map
tmap_mode("plot")
tmap_arrange(districts_choropleth)

More specifically, we can observe that conflict hotspots are mainly found in the districts of Yinmarbin, Shwebo, Pakokku and Mandalay which lies in the central regions of Myanmar.

4.2 Identifying States with Highest Proportion of Conflicts

Instead, let us also explore the top 10 states with the highest proportions of armed conflict events.

Count number of conflicts by states
conflict_count <- conflict_data_sf %>%
  group_by(ST) %>%
  summarise(total_count_ST = n()) %>%
  st_drop_geometry() %>%
  select(ST, total_count_ST)

# Perform the join
boundary_sf <- boundary_sf %>%
  left_join(conflict_count, by = "ST")

Likewise, I’ll add a new column called proportion_ST to represent the proportion based on each Myanmar state.

boundary_sf <- boundary_sf %>%
  mutate(proportion_ST = total_count_ST / sum(total_count_ST))

head(boundary_sf[c('ST','total_count_ST','proportion_ST')])
Simple feature collection with 6 features and 3 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -14915.04 ymin: 1736124 xmax: 187961.7 ymax: 2051144
Projected CRS: WGS 84 / UTM zone 47N
          ST total_count_ST proportion_ST                       geometry
1 Ayeyarwady            853            NA MULTIPOLYGON (((90859.89 20...
2 Ayeyarwady            853            NA MULTIPOLYGON (((75991.51 17...
3 Ayeyarwady            853            NA MULTIPOLYGON (((115559 1928...
4 Ayeyarwady            853            NA MULTIPOLYGON (((39919.39 18...
5 Ayeyarwady            853            NA MULTIPOLYGON (((-6302.348 1...
6 Ayeyarwady            853            NA MULTIPOLYGON (((93411.72 17...

At a quick glance, we can see that central and southern parts of Myanmar have the highest proportions of armed conflict events occurring, namely in Sagaing, Mandalay, Magway and Yangon states as indicated in the map below (darkest shade of blue).

Create the points map
states_choropleth <-
tm_shape(boundary_sf) +
  tm_fill("proportion_ST",
          n=10,
          title="Proportion",
          style="equal",
          palette="Blues") +
  tm_borders(lwd=0.2,
             alpha=1) +
  tm_text(text = "ST", 
          size = 0.2, 
          col = "black",
          fontface = "bold") +
  tm_layout(main.title = "Distribution of Conflict Points Across States",
            legend.outside=FALSE,
            main.title.size=1)
# Plot the map
tmap_mode("plot")
tmap_arrange(states_choropleth)

4.3 Standard Distances of Top 4 Districts

Next, I’ve plot the standard distances for the four districts that I’m most interested in. The codes here will calculate the average location of the conflict points within the district (i.e. the mean centre) and measure the dispersion of points around the mean centre.

The circle drawn around the mean centre shows us the standard distance to which we can observe that Mandalay has the smallest standard distance and hence, has the most conflict points closely clustered around the mean centre.

Find hot spot of Yinmarbin District
# Define the districts
districts <- c("Yinmarbin", "Shwebo", "Pakokku", "Mandalay")

for (district in districts) {
  # Filter the conflict data for the current district
  conflict_district <- filter(conflict_data_sf, DT == district)
  boundary_district <- filter(boundary_sf, DT == district)

  # Create a combined ppp and owin object
  district_owin <- as.owin(boundary_district)
  ppp_obj <- as.ppp(st_geometry(conflict_district))

  # Handle duplicates
  ppp_obj <- rjitter(ppp_obj, retry = TRUE, nsim = 1, drop = TRUE)

  # Mask ppp object with the boundary
  district_ppp_owin <- ppp_obj[district_owin]
  district_ppp_owin <- rescale(district_ppp_owin, 1000, "km")

  # Calculate the mean center
  mean_coords <- c(mean(district_ppp_owin$x), mean(district_ppp_owin$y))

  # Calculate the standard distance
  sd_x <- sd(district_ppp_owin$x)
  sd_y <- sd(district_ppp_owin$y)
  standard_distance <- sqrt(sd_x^2 + sd_y^2)

  # Prepare to plot for the current district
  plot(district_ppp_owin, 
       main = paste("Standard Distance in", district), 
       pch = 19, 
       col = "lightgrey",
       xlim = range(district_ppp_owin$x) + c(-5, 5),  # Expand limits
       ylim = range(district_ppp_owin$y) + c(-5, 5))  # Expand limits

  # Highlight the mean center with an 'X'
  points(mean_coords[1], mean_coords[2], pch = 4, col = "red", cex = 1)  # Smaller 'X'

  # Create a circle around the mean center for the standard distance
  bearing <- seq(0, 2 * pi, length.out = 360)
  circle_x <- mean_coords[1] + standard_distance * cos(bearing)
  circle_y <- mean_coords[2] + standard_distance * sin(bearing)

  # Draw the circle
  lines(circle_x, circle_y, col = 'red', lwd = 2)
}

# Reset the plot layout
par(mfrow = c(1,1), mar = c(0,0,1,0)) 

5. 1st Order Spatial Point Patterns Analysis

5.1 Kernel Density Estimation

5.1.1 Working with Fixed Bandwidth Methods

Using the geospatial data sets prepared, I will now perform 1st order spatial point pattern analysis by leveraging kernel density estimation (KDE) to understand the intensity of conflicts in different regions.

I will be using a variety of fixed bandwidth methods via density() of the spatstat package, to determine the most optimal method for this analysis. Namely using bw.diggle(), bw.ppl(), bw.CvL() and bw.scott().

Steps taken to calculate the KDE:

  1. Extract the masked ppp object for the current quarter.
  2. Compute the kernel density estimate by setting the signma parameters.
  3. Plot the kernel density estimate using plot() where “Bw” represents the optimal bandwidth

For the purposes of identifying the most optimal bandwidth method, I will create a ppp_obj using the 2021 Q1 conflict events first to assist my decision-making.

# Set Up
ppp_obj = masked_ppp_list_km$`2021 Q1`
colours <- colorRampPalette(c("midnightblue", "skyblue"))(100)

The bw.diggle() bandwidth is referred to as Diggle’s cross-validation bandwidth which minimises the mean-squared error (MSE) to balance between under and over-smoothing.

# bw.diggle()
kde_conflict_bw_diggle <- density(ppp_obj,
                           sigma=bw.diggle,
                           edge=TRUE,
                           kernel="gaussian")
optimal_bw_d = floor(bw.diggle(ppp_obj)[[1]]*10)/10
plot(kde_conflict_bw_diggle, main = paste("BW: diggle", "(",optimal_bw_d,"km)"), col = colours)

The second bandwidth method I attempted using is bw.ppl(), This method chooses the bandwidth that minimises the likelihood cross-validation score and improving the prediction accuracy of the kernel density estimate.

# bw.ppl()
kde_conflict_bw_ppl <- density(ppp_obj,
                           sigma=bw.ppl,
                           edge=TRUE,
                           kernel="gaussian")
optimal_bw_p = floor(bw.ppl(ppp_obj)[[1]]*10)/10
plot(kde_conflict_bw_ppl, main = paste("Bw: ppl", "(",optimal_bw_p,"km)"), col = colours)

Thirdly, let’s explore the bandwidth method bw.CvL(), also known as Cronie and Van Lieshout cross-validation, designed to provide an optimal, adaptive bandwidth for inhomogeneous point patterns.. Similar to bw.ppl(), it aims to reduce the error measure but also aims to balance over and under-fitting based on the spatial structure of the data.

# bw.CvL()
kde_conflict_bw_CvL <- density(ppp_obj,
                           sigma=bw.CvL,
                           edge=TRUE,
                           kernel="gaussian")
optimal_bw_c = floor(bw.CvL(ppp_obj)[[1]]*10)/10
plot(kde_conflict_bw_CvL, main = paste("Bw: CvL (",optimal_bw_c,"km)"), col = colours)

Lastly, I will explore the bw.scott() bandwidth method. This method returns separate bandwidths for the x- and y-axes which is ideal for our spatial data that contains both x and y components. I will combine these bandwidths into a single value for isotropic kernel density estimation by taking the taking the geometric mean as shown in the value returned by sigma_combined.

# bw.scott()
bw_values <- bw.scott(ppp_obj)
sigma_x <- bw_values[1]
sigma_y <- bw_values[2]
sigma_combined <- sqrt(sigma_x * sigma_y)

kde_conflict_bw_scott <- density(ppp_obj,
                           sigma = sigma_combined,
                           edge = TRUE,
                           kernel = "gaussian")

optimal_bw_s = floor(sigma_combined*10)/10
plot(kde_conflict_bw_scott, main = paste("Bw: scott", "(",optimal_bw_s,"km)"), col = colours)

Selecting a Bandwidth Method

Pramanik N. (2019) proposes that choosing an optimal bandwidth method is crucial in fitting the data appropriately to balance between bias and variance. As such, the four methods listed below cater to different types of data depending on how varied the densities are spread across and the granularity of conflict events. Let’s explore which bandwidth method is optimal for our dataset.

  • bw.diggle(): appears to be effective for homogeneous data in seeing general conflict hotspots and uses the smallest bandwidth size
  • bw.ppl(): tends to choose slightly larger bandwidths, it provide more localised density estimates which highlights finer spatial details. As such, we can see more variability and finer details in the density distribution, with more variation between high- and low-density areas.
  • bw.CvL(): the kernel density plot shows that CvL makes a good attempt in balancing between detail and smoothness, making it more suitable for capturing the overall density trends in spatial data with some local structures highlighted.
  • bw.scott(): as shown, the geometric mean ensures equal smoothing in both x and y directions, and is largely similar to bw.Cvl(), making it a good choice for a balanced and general and fast overview of the spatial data distribution.
Plot all bandwidth methods
par(mfrow = c(2,2), mar = c(0,0,1,0)) 
plot(kde_conflict_bw_diggle, main = paste("diggle (",optimal_bw_d,"km)"), col = colours)
plot(kde_conflict_bw_ppl, main = paste("ppl (",optimal_bw_p,"km)"), col = colours)
plot(kde_conflict_bw_CvL, main = paste("CvL (",optimal_bw_c,"km)"), col = colours)
plot(kde_conflict_bw_scott, main = paste("scott (",optimal_bw_s,"km)"), col = colours)

# Plot histogram to compare KDE
par(mar = c(2,2,2,2),mfrow = c(2,2))
hist(kde_conflict_bw_diggle, main = "diggle")
hist(kde_conflict_bw_ppl, main = "ppl")
hist(kde_conflict_bw_CvL, main = "scott")
hist(kde_conflict_bw_scott, main = "CvL")

💡 Decision: I decided to use bw.CvL() for computing the KDE of the masked ppp objects based on each quarter. The KDE function outputs a relatively smooth density estimate that isn’t too detailed like bw.ppl() and not as generalised as bw.scott().

This coincides with wider range of spatial concentrations captured by the histogram, doing a good job with both capturing the non-homogeneous data spread of our Myanmar conflict data. Moreover, it isn’t as computationally heavy as bw.ppl().

Putting Together our Fixed KDE using bw.CvL()

Now, let us perform the KDE computation for the conflict events across all quarters using bw.CvL().

# Calculate density using bw.CvL()
par(mfrow = c(2,3), mar = c(0,0,1,0)) 
for (quarter in names(masked_ppp_list_km)) {
  ppp_obj = masked_ppp_list_km[[quarter]]
  kde_conflict_bw <- density(ppp_obj,
                             sigma=bw.CvL,
                             edge=TRUE,
                             kernel="gaussian")
  optimal_bw = floor(bw.CvL(ppp_obj)[[1]]*10)/10
  plot(kde_conflict_bw, main = paste(quarter, "(Bw:",optimal_bw,"m)"), col = colours)
}

Observations

Bandwidth Size: A bandwidth of around 60,000 to 100,000 is considered relatively large as compared to the bandwidth returned from using bw.diggle() and bw.ppl(). Hence, this results in a smoother density estimate with less emphasis on local clusters as indicated in the generalised spatial trends.

Note: Using a fixed KDE is beneficial for our quarterly analysis but the bandwidth returned is different for each quarter. Later on, I’ll take the average bandwidth across different time periods so we can focus on spatial distribution over time.

5.1.2 Working with Different Kernel Methods

Bias is an inevitable issue when it comes to choosing the type of kernel function. Rajagopalan B., et al (1997) states that kernels such as Epanechnikov, Quartic and Disc suffer a “leakage problem” where some of the probability mass is truncated at the boundary, leading to boundary bias. At the same time, the Gaussian kernel has no truncation, but small “leakage” beyond data range of the boundary inadvertantly causes bias to exist still. Shen et al. (2020) further proposes that the bandwidth is a more of an influential factor in kernel estimation than the selection of different kernel functions.

Let us evaluate if this is true by experimenting with a variety of kernel methods as they ultimately control how we weight points within the bandwidth radius. The default kernel in density.ppp() is the gaussian. alternatives such as epanechnikov, quartic, and disc are also available. I will use the 2021 Q1 conflict data to assist in identifying the most optimal kernel method.

# Set Up
par(mfrow = c(2,2), mar = c(0,0,1,0)) 
ppp_obj = masked_ppp_list_km$`2021 Q1`

# Using the gaussian kernel
kde_conflict_g <- density(ppp_obj,
                           sigma=bw.CvL,
                           edge=TRUE,
                           kernel="gaussian")
plot(kde_conflict_g, main="Gaussian Kernel", col = colours)

# Using the epanechniko kernel
kde_conflict_e <- density(ppp_obj,
                           sigma=bw.CvL,
                           edge=TRUE,
                           kernel="epanechnikov")
plot(kde_conflict_e, main="Epanechnikov Kernel", col = colours)

# Using the quartic kernel
kde_conflict_q <- density(ppp_obj,
                           sigma=bw.CvL,
                           edge=TRUE,
                           kernel="quartic")
plot(kde_conflict_q, main="Quartic Kernel", col = colours)

# Using the disc kernel
kde_conflict_d <- density(ppp_obj,
                           sigma=bw.CvL,
                           edge=TRUE,
                           kernel="disc")
plot(kde_conflict_d, main="Disc Kernel", col = colours)

# Plot histogram to compare KDE
par(mar = c(2,2,2,2),mfrow = c(2,2))
hist(kde_conflict_g, main = "Gaussian")
hist(kde_conflict_e, main = "Epanechnikov")
hist(kde_conflict_q, main = "Quartic")
hist(kde_conflict_d, main = "Disc")

Observations

As expected, we don’t see major variations in the smoothness and spread of KDE across a range of distances, and bias continues to be present across all kernel methods as seen from the long-tails. However, there are slight differences in the sharpness of the bandwidth radius in terms of how localised or widespread our data points are being captured.

  • Gaussian: provides a localised density estimate over the entire spatial extent as compared to epanechnikov and quartic. It is good at highlighting variance and opposing ends of conflict intensities as shown by the wider range used in the legend.
  • Epanechnikov: It is more efficient than the gaussian in terms of variance but produces a slightly rougher surface. It is also more localised than the quartic kernel, focusing on areas near each point, with a sharper boundary at the bandwidth limit.
  • Quartic: Results in a good balance between smoothness and localised influence, smoother than epanechnikov but with similar properties. It appears suitable for moderate smoothing and sharper focus on local patterns.
  • Disc: results in the sharpest density estimate as compared to the other three kernels as all points within a certain distance are made to have equal influence and zero influence beyond that distance.

Decision: Hence, I will use the quartic kernel method to ensure a relatively smooth density estimate with emphasis on local points over distant ones.

5.1.3.1 Using kernel = ‘quartic’

As such, I run the density estimate computation using kernel = 'quartic'.

# Using 'quartic' kernel
par(mfrow = c(2,3), mar = c(0,0,1,0)) 

for (quarter in names(masked_ppp_list_km)) {
  ppp_obj = masked_ppp_list_km[[quarter]]
  kde_conflict_bw <- density(ppp_obj,
                             sigma=bw.CvL,
                             edge=TRUE,
                             kernel="quartic")
  optimal_bw = floor(bw.CvL(ppp_obj)[[1]]*1000)/1000
  plot(kde_conflict_bw, main = paste(quarter, "(Bw:",optimal_bw,"km)"), col = colours)
}

Observations

We can see high densities of armed conflict in the central and southern regions of Myanmar but more can be uncovered from conflict data. Let’s proceed to the next section.

5.1.3.2 Using sigma = 71.831

For all subsequent fixed KDE computations, I will assign sigma using the average of the CvL bandwidth returned from each quarter. Here’s the calculations based on the plots returned:

Average bandwidth size = 61.649 + 64.386 + 74.501 + 114.08 + 103.863 + 103.863 + 95.567 + 56.757 + 57.752 + 66.968 + 49.649 + 48.135 + 48.135 + 60.323) / 14 = 71.831

Let’s recompute the Fixed KDE based on the newly calculated average bandwidth such that sigma = 71.831. I’ll store the quarterly KDE outputs into a list called kde_conflict_bw_list.

# Add KDE into this list
kde_conflict_bw_list <- list()
for (quarter in names(masked_ppp_list_km)) {
  ppp_obj = masked_ppp_list_km[[quarter]]
  kde_conflict_bw <- density(ppp_obj,
                             sigma=71.831,
                             edge=TRUE,
                             kernel="quartic")
  kde_conflict_bw_list[[quarter]] <- kde_conflict_bw
}
# Plot graph
par(mfrow = c(2,3), mar = c(0,0,1,0)) 

for (quarter in names(kde_conflict_bw_list)){
  kde_conflict_bw <- kde_conflict_bw_list[quarter]
  plot(kde_conflict_bw, main = paste(quarter, "(Bw: 71.831 km)"), col = colours)
}

5.1.3 Working with Adaptive KDE

As seen above, fixed bandwidths tend to oversmooth the mode of the distribution. On the contrary, the adaptive kernel estimate has the ability to reduce variability of estimates in areas with low density and increases it in areas with higher density (The Stata Journal, 2003).

Once again, let us use the 2021 Q1 conflict data to illustrate the difference in outputs of all three adaptive methods.

# Set Up
ppp_obj = masked_ppp_list_km$`2021 Q1`

# Using Voronoi adaptive KDE
vd_adaptive_kde <- adaptive.density(ppp_obj, f=1, method="voronoi")

# Plot
par(mar = c(0,1,1,1))
plot(vd_adaptive_kde, main = "Voronoi-Dirichlet Adaptive KDE", col = colours)

# Using adaptive KDE
adaptive_kde <- adaptive.density(ppp_obj, method="kernel")

# Plot
par(mar = c(0,1,1,1))
plot(adaptive_kde, main = "Adaptive KDE", col = colours)

I’ve used a relatively larger number of neighbours (i.e. k = 10) to provide a smoother, more general density estimate to capture broader trends and may smooth out details.

# Using nearest neighbour adaptive KDE
nn_kde <- nndensity(ppp_obj, k=10)

# Plot
par(mar = c(0,1,1,1))
plot(nn_kde, main = "Nearest Neighbour Adaptive KDE", col = colours)

We can also compare the performance of each method based on the top 4 states with highest proportion of conflicts as highlighted earlier.

Comparing the three Adaptive KDE Types

From the outputs above, it appears that there is no major differences between the distribution of KDE values returned across the three methods, where there is high concentration of points in a specific area. Hence, we will choose to go with Adapative Kernel method.

Let’s compare the results of my two selected fixed and adaptive KDEs (E.g. Magway District)

We can observe how adaptive kernels provides a more detailed picture of conflict spatial distribution but since it’s largely localised, conflict spots require more effort in identifying and can be computationally heavy for this exercise.

Additionally, varying bandwidth makes comparisons across regions or time periods (like quarters) more difficult because the scale of smoothing is not constant across space and time.

5.1.4 Converting Gridded KDE Output into Raster

Next, we need to convert the KDE output to KDE raster layers before it can be viewed using tmap.

Step 1) Converting KDE to Spatial Grid Data Frame

library(grid)
plot_list <- list()
for (quarter in names(kde_conflict_bw_list)) {
  ppp_obj <- kde_conflict_bw_list[[quarter]]
  gridded_ppp_obj <- as(ppp_obj, "SpatialGridDataFrame")
  plot_list[[quarter]] <- spplot(gridded_ppp_obj, main = paste(quarter), col.regions = colours)
}

library(gridExtra)
plot_list_subset1 <- plot_list[1:6]
plot_list_subset2 <- plot_list[7:12]
plot_list_subset3 <- plot_list[13:14]
grid.newpage()
grid.arrange(grobs = plot_list_subset1, ncol = 3, nrow = 2)
grid.newpage() 
grid.arrange(grobs = plot_list_subset2, ncol = 3, nrow = 2)
grid.newpage()
grid.arrange(grobs = plot_list_subset3, ncol = 3, nrow = 2)

Step 2) Rasterisation of Grid Outputs & Assigning Projection Systems

gridded_ppp_obj_raster_list <- list()
for (quarter in names(kde_conflict_bw_list)) {
  gridded_ppp_obj = kde_conflict_bw_list[[quarter]]
  gridded_ppp_obj_raster <- raster(gridded_ppp_obj)
  projection(gridded_ppp_obj_raster) <- CRS("+init=EPSG:32647")
  gridded_ppp_obj_raster_list[[quarter]] <- gridded_ppp_obj_raster
}

# Inspect for 2024 Q2
gridded_ppp_obj_raster_list$`2024 Q2`
class      : RasterLayer 
dimensions : 128, 128, 16384  (nrow, ncol, ncell)
resolution : 7.302001, 16.30032  (x, y)
extent     : -210.0086, 724.6476, 1072.026, 3158.467  (xmin, xmax, ymin, ymax)
crs        : +proj=utm +zone=47 +datum=WGS84 +units=m +no_defs 
source     : memory
names      : layer 
values     : 7.061861e-06, 0.01847887  (min, max)

Step 3) Plot Maps

tmap_mode("plot")
plots_by_quarter <- list()
for (quarter in names(gridded_ppp_obj_raster_list)){
  gridded_ppp_obj_raster = gridded_ppp_obj_raster_list[[quarter]]
  raster_plot <- tm_shape(gridded_ppp_obj_raster) +
    tm_raster("layer", title="Density", palette="Blues") +
    tm_layout(legend.position = c("left","bottom"),frame=FALSE, main.title = quarter,
            main.title.size=1, main.title.position = "center", legend.text.size = 0.5,legend.title.size = 0.7) 
  plots_by_quarter[[quarter]] <- raster_plot
}

tmap_arrange(plots_by_quarter[1:14], ncol=5, nrow=3)

Observations

Plotting raster grid versions of KDE outputs uses discrete colour ranges which does a good job in highlighting gradual changes in conflict events across an area. Since 2021 Q2, more conflicts are seen in Southern parts of Myanmar. The density range differs for each quarter but we can see an increase in no. of armed conflicts per kilometre from 2021 Q1 to 2022 Q2, which stagnates in density and increases again in 2023 Q4.

5.2 Nearest Neighbour Analysis

Our current analyses does not reveal patterns of clustering or dispersion, to which Michael J. Crawley proposes to employ Clark-Evans test spatial randomness for its simplicity and applicability for first-order spatial analysis, which means checking for overall spatial randomness based on nearest-neighbor distances. (Crawley M. J. , 2007)

Clark-Evans Test

The test checks whether the observed point pattern of armed conflicts in Myanmar shows clustering (points are closer than expected under randomness), dispersion (points are more spread out), or randomness.

5.2.1 Clark-Evans Test (Myanmar)

The test hypotheses are:

  • Ho = The distribution of armed conflicts in Myanmar are randomly distributed.
  • H1= The distribution of armed conflicts in Myanmar are not randomly distributed.
  • The 95% confident interval will be used.

We will conduct the test using clarkevans.test() of statspat.

for (quarter in names(masked_ppp_list_km)) {
  ppp_obj = masked_ppp_list_km[[quarter]]
  print(quarter)
  print(clarkevans.test(ppp_obj,
                  correction="none",
                  clipregion="boundary_sf",
                  alternative=c("clustered"),
                  nsim=99))
}
[1] "2024 Q2"

    Clark-Evans test
    No edge correction
    Z-test

data:  ppp_obj
R = 0.26663, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)

[1] "2024 Q1"

    Clark-Evans test
    No edge correction
    Z-test

data:  ppp_obj
R = 0.23563, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)

[1] "2023 Q4"

    Clark-Evans test
    No edge correction
    Z-test

data:  ppp_obj
R = 0.21795, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)

[1] "2023 Q3"

    Clark-Evans test
    No edge correction
    Z-test

data:  ppp_obj
R = 0.22002, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)

[1] "2023 Q2"

    Clark-Evans test
    No edge correction
    Z-test

data:  ppp_obj
R = 0.24485, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)

[1] "2023 Q1"

    Clark-Evans test
    No edge correction
    Z-test

data:  ppp_obj
R = 0.24365, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)

[1] "2022 Q4"

    Clark-Evans test
    No edge correction
    Z-test

data:  ppp_obj
R = 0.22139, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)

[1] "2022 Q3"

    Clark-Evans test
    No edge correction
    Z-test

data:  ppp_obj
R = 0.23974, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)

[1] "2022 Q2"

    Clark-Evans test
    No edge correction
    Z-test

data:  ppp_obj
R = 0.22989, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)

[1] "2022 Q1"

    Clark-Evans test
    No edge correction
    Z-test

data:  ppp_obj
R = 0.21976, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)

[1] "2021 Q4"

    Clark-Evans test
    No edge correction
    Z-test

data:  ppp_obj
R = 0.21341, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)

[1] "2021 Q3"

    Clark-Evans test
    No edge correction
    Z-test

data:  ppp_obj
R = 0.21808, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)

[1] "2021 Q2"

    Clark-Evans test
    No edge correction
    Z-test

data:  ppp_obj
R = 0.17458, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)

[1] "2021 Q1"

    Clark-Evans test
    No edge correction
    Z-test

data:  ppp_obj
R = 0.24696, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)
Observations

For a 95% confidence level, If the p-value < 0.05, I will reject the null hypothesis of complete spatial randomness and check if data is uniform (R > 1) or clustered (R < 1).

With that said, all tests conducted across each quarter rejects the null hypothesis as p < 0.05 and spatial points are found to be clustered since R < 1.

5.2.2 Clark-Evans Test (Top 4 Districts)

# Yinmarbin District
clarkevans.test(yinmarbin_ppp_owin,
                  correction="none",
                  clipregion="boundary_yinmarbin",
                  alternative=c("clustered"),
                  nsim=39)

    Clark-Evans test
    No edge correction
    Z-test

data:  yinmarbin_ppp_owin
R = 0.13718, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)
# Shwebo District
clarkevans.test(shwebo_ppp_owin,
                  correction="none",
                  clipregion="boundary_shwebo",
                  alternative=c("clustered"),
                  nsim=39)

    Clark-Evans test
    No edge correction
    Z-test

data:  shwebo_ppp_owin
R = 0.20141, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)
# Pakokku District
clarkevans.test(pakokku_ppp_owin,
                  correction="none",
                  clipregion="boundary_pakokku",
                  alternative=c("clustered"),
                  nsim=39)

    Clark-Evans test
    No edge correction
    Z-test

data:  pakokku_ppp_owin
R = 0.20698, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)
# Mandalay District
clarkevans.test(mandalay_ppp_owin,
                  correction="none",
                  clipregion="boundary_mandalay",
                  alternative=c("clustered"),
                  nsim=39)

    Clark-Evans test
    No edge correction
    Z-test

data:  mandalay_ppp_owin
R = 0.16014, p-value < 2.2e-16
alternative hypothesis: clustered (R < 1)
Observations

I will reject the null hypothesis of complete spatial randomness since the p values of each district is smaller than 0.05. Additionally, we can observe clustering of spatial points since R < 1 is returned for all districts. However, this test isn’t sufficient in highlighting the statistical significance of the spatial patterns - I’ll use the Monte Carlo Simulation later to handle this.

5.3 Further Data Exploration

By using the fixed KDE with CvL bandwidth and quartic kernel, let’s see what insights can we glean from the density of conflicts in Myanmar.

5.3.1 KDE by Event Type

First, let’s identify the unique event types in this dataset.

# Check unique events
unique(conflict_data_sf$event_type)
[1] "Battles"                    "Strategic developments"    
[3] "Violence against civilians" "Explosions/Remote violence"

Now, let us analyse the kernel density estimate of each unique event type found in conflict_data_sf to identify hot and cold spots across Myanmar.

Plot the KDE based on Event Type
# Set Up
par(mfrow = c(2,2), mar = c(0,0,1,0)) 

conflict_data_sf %>%
  group_by(event_type) %>%
  group_split() -> conflict_by_event_type

# Convert the sf object to owin
district_boundary <- as.owin(st_as_sfc(boundary_sf))

kde_list <- lapply(seq_along(conflict_by_event_type), function(i) {
  data <- conflict_by_event_type[[i]]
  event_type <- unique(data$event_type)
  ppp_obj <- as.ppp(st_geometry(data), W = district_boundary)
  ppp_obj <- rescale(ppp_obj, 1000, "km")
  kde <- density(ppp_obj,
                 sigma=71.831,
                 edge=TRUE,
                 kernel="quartic")
  plot(kde, main = paste(event_type), col=colours)
  return(kde)
})

Observations

We can almost see an equal spread of all four event types, with explosions and strategic violence being more dominantly found in Central-Western Myanmar, in the Sagaing and Mandalay states, followed by battles and violence against citizens. All events seem to plaque the Southern states (e.g. Yangon) with the exception of battles.

5.3.2 KDE Across Top 4 States With Most Conflicts

Previously, we identified the top 4 states with the highest proportions of conflicts as Sagaing, Mandalay, Magway and Yangon. We can delve deeper into each state by analysing the intensity of conflicts across these states using density().

Plot the KDE of Top 4 States
# Set Up
par(mfrow = c(2,2), mar = c(0,0,1,0))

# Sagaing
district_boundary <- as.owin(st_as_sfc(boundary_sagaing))
ppp_obj <- as.ppp(st_geometry(conflict_sagaing), W = district_boundary)
ppp_obj_sagaing <- rescale(ppp_obj, 1000, "km")
kde <- density(ppp_obj_sagaing,
               sigma=71.831,
               edge=TRUE,
               kernel="quartic")
plot(kde, main = paste("Sagaing"), col=colours)

# Mandalay
district_boundary <- as.owin(st_as_sfc(boundary_mandalay))
ppp_obj <- as.ppp(st_geometry(conflict_mandalay), W = district_boundary)
ppp_obj_mandalay <- rescale(ppp_obj, 1000, "km")
kde <- density(ppp_obj_mandalay,
               sigma=71.831,
               edge=TRUE,
               kernel="quartic")
plot(kde, main = paste("Mandalay"), col=colours)

# Magway
district_boundary <- as.owin(st_as_sfc(boundary_magway))
ppp_obj <- as.ppp(st_geometry(conflict_magway), W = district_boundary)
ppp_obj_magway <- rescale(ppp_obj, 1000, "km")
kde <- density(ppp_obj_magway,
               sigma=71.831,
               edge=TRUE,
               kernel="quartic")
plot(kde, main = paste("Magway"), col=colours)

# Yangon
district_boundary <- as.owin(st_as_sfc(boundary_yangon))
ppp_obj <- as.ppp(st_geometry(conflict_yangon), W = district_boundary)
ppp_obj_yangon <- rescale(ppp_obj, 1000, "km")
kde <- density(ppp_obj_yangon,
               sigma=71.831,
               edge=TRUE,
               kernel="quartic")
plot(kde, main = paste("Yangon"), col=colours)

Observations

It’s interesting that armed conflict isn’t evenly distributed across the states though it does seem that armed conflict has inflicted the entire state of Yangon. Nonetheless, it is worth noting that Yangon is relatively smaller in size than the other three states and that will increase the density of conflict quite significantly.

5.3.3 KDE of Top 4 States by Event Type

It’ll also be interesting to breakdown each top 4 state by the event type category as shown.

Observations

The kernel density of violence against civillians is generally found to be the lowest amongst all conflict events. Additionally, all types of armed conflicts tend to occur repeatedly in the same parts of each state. E.g. conflicts regarding strategic development tend to happen in Southern part of the Sagaing state, just as it is for explosions/remote violence.

5.5.4 KDE by Interaction Type

# Convert Interaction Code to Text
unique(conflict_data_sf$interaction)
 [1] 13 10 12 37 70 17 22 11 33 30 27 20 23 60 47 28 38 80 78 18 24 14
library(dplyr)

# Create a named vector for mapping interaction values
interaction_map <- c(
  "10" = "SOLE STATE FORCES ACTION",
  "11" = "STATE FORCES VERSUS STATE FORCES",
  "12" = "STATE FORCES VERSUS REBELS",
  "13" = "STATE FORCES VERSUS POLITICAL MILITIA",
  "14" = "STATE FORCES VERSUS IDENTITY MILITIA",
  "17" = "STATE FORCES VERSUS CIVILIANS",
  "18" = "STATE FORCES VERSUS EXTERNAL/OTHER FORCES",
  "20" = "SOLE REBEL ACTION",
  "22" = "REBELS VERSUS REBELS",
  "23" = "REBELS VERSUS POLITICAL MILITIA",
  "24" = "REBELS VERSUS IDENTITY MILITIA",
  "27" = "REBELS VERSUS CIVILIANS",
  "28" = "REBELS VERSUS OTHERS",
  "30" = "SOLE POLITICAL MILITIA ACTION",
  "33" = "POLITICAL MILITIA VERSUS POLITICAL MILITIA",
  "34" = "POLITICAL MILITIA VERSUS IDENTITY MILITIA",
  "37" = "POLITICAL MILITIA VERSUS CIVILIANS",
  "38" = "POLITICAL MILITIA VERSUS OTHERS",
  "40" = "SOLE IDENTITY MILITIA ACTION",
  "44" = "IDENTITY MILITIA VERSUS IDENTITY MILITIA",
  "47" = "IDENTITY MILITIA VERSUS CIVILIANS",
  "48" = "IDENTITY MILITIA VERSUS OTHER",
  "70" = "SOLE CIVILIAN ACTION",
  "77" = "CIVILIANS VERSUS CIVILIANS",
  "78" = "OTHER ACTOR VERSUS CIVILIANS",
  "80" = "SOLE OTHER ACTION",
  "88" = "OTHER VERSUS OTHER"
)

conflict_data_sf <- conflict_data_sf %>%
  rename(interaction_code = interaction) %>%
  mutate(interaction = interaction_map[as.character(interaction_code)])  

# View the updated dataframe
head(conflict_data_sf[c("interaction_code","interaction")])
Simple feature collection with 6 features and 2 fields
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: 76997.72 ymin: 2428487 xmax: 214961 ymax: 2533434
Projected CRS: WGS 84 / UTM zone 47N
# A tibble: 6 × 3
  interaction_code interaction                                   geometry
             <dbl> <chr>                                      <POINT [m]>
1               13 STATE FORCES VERSUS POLITICAL MILI…   (214961 2452068)
2               13 STATE FORCES VERSUS POLITICAL MILI… (198303.2 2499463)
3               13 STATE FORCES VERSUS POLITICAL MILI… (189105.4 2533434)
4               10 SOLE STATE FORCES ACTION            (160913.9 2522331)
5               13 STATE FORCES VERSUS POLITICAL MILI…   (146213 2428487)
6               10 SOLE STATE FORCES ACTION            (76997.72 2447719)
1. Involvement of Political Militia

The conflict landscape in Myanmar is complex and fragmented, with political militias focused on advancing political agendas or controlling governance structures. Let’s discover how the density of armed conflicts differ for each interaction of the political militia. They can be associated as either of these groups.

  • Pro-Government Militias: supporting the central government or Tatmadaw, such as Pyithu Sit and BGFs.
  • Pro-Democracy Militias: emerged as part of the resistance movement against the military junta. E.g. the People’s Defence Forces (PDF)
Plot conflicts involving political militia
codes <- c('13','23','33')
par(mfrow = c(3,1), mar = c(0,0,1,0)) 

for(code in codes) {
  # Filter for the specific interaction type
  filtered_data <- conflict_data_sf %>%
    filter(interaction_code == code)
  
  # Convert to ppp object
  boundary_owin <- as.owin(boundary_sf)
  ppp_owin_obj <- as.ppp(filtered_data$geometry, W = boundary_owin)
  ppp_owin_obj_km <- rescale(ppp_owin_obj, 1000, "km")
  
  # Compute KDE
  kde <- density(ppp_owin_obj_km,
                 sigma=71.831,
                 edge=TRUE,
                 kernel="quartic")
  
  # Plot the KDE
  plot(kde, main = paste(interaction_map[code]), col = colours)
}

Observations

The political army or organisations tend to be less involved with unarmed civilians but state forces, rebels and other policial militia.

  1. Against State Forces
  • High density of conflicts with state forces (e.g. police) in central and western Myanmar
  • These regions have significant ethnic diversity and political grievances. E.g. Rakhine State in the west is home to ethnic Rakhine groups who have historically sought greater autonomy and have been involved in conflicts with state forces.
  • Central regions, such as Sagaing and Magway, have also seen intense resistance against state forces, particularly after the 2021 military coup.
  1. Against Rebels
  • Rebels and political militias operate throughout Myanmar with varying objectives but of a much smaller intensity than with other actors (seen from its low density scale).
  • Clashes between rebels and political militias can stem from differing political goals, territorial disputes, or strategic interests.
  • For instance, a rebel group might view a political militia as a rival in a contested region or as a competitor for local support.
  1. Against other Political Miliitia
  • In central and western Myanmar, where political militias are more active, there is often competition among different resistance groups. These militias may have different leadership, ideologies, or goals, leading to clashes.
2. Involvement of Identity Militia

Identity militias are typically organised around ethnic, religious, or cultural identities. heir main goal is to protect the interests, rights, and autonomy of specific identity-based groups within Myanmar.

Plot conflicts involving identity militia
codes <- c('14','24')
par(mfrow = c(2,1), mar = c(0,0,1,0)) 

for(code in codes) {
  # Filter for the specific interaction type
  filtered_data <- conflict_data_sf %>%
    filter(interaction_code == code)
  
  # Convert to ppp object
  boundary_owin <- as.owin(boundary_sf)
  ppp_owin_obj <- as.ppp(filtered_data$geometry, W = boundary_owin)
  ppp_owin_obj_km <- rescale(ppp_owin_obj, 1000, "km")
  
  # Compute KDE
  kde <- density(ppp_owin_obj_km,
                 sigma=71.831,
                 edge=TRUE,
                 kernel="quartic")
  
  # Plot the KDE
  plot(kde, main = paste(interaction_map[code]), col = colours)
}

Observations

Armed conflicts involving the identity militia is of a much smaller intensity than political militia with about 5e^-0.5 or 3 conflict events per kilometre.

  1. State Forces vs Identity Militia
  • Conflicts between these two parties are concerntrated in the Western region which has a history of conflict, including uprisings against colonial and post-colonial governments.

  • Additionally, the region’s proximity to Bangladesh has some influence on local conflicts as ethnic groups in Myanmar may have connections or support from across the border.

  • In the Rakhine State of western Myanmar, historical grievances and demands for autonomy has also fueled tensions between these groups and the central government (Crisis Group, 2024).

  1. Rebels vs identity militia
  • Higher intensity of armed conflicts are seen in North Myanmar
  • Northern Myanmar hosts various ethnic armed groups and rebel factions that seek autonomy or independence. (E.g. Kachin, Shan and Chin)
  • Northern Myanmar is rich in natural resources, e.g. jade and timber. Control over these resources can be a significant driver of conflict, as various armed groups vie for control and when there’s a lack of economic opportunities. (Fishbein & Lusan, 2022)
  • Unresolved conflicts from colonial and post-colonial periods continue to influence current tensions between different ethnic groups and armed factions in the Northen region.
3. Involvement of Civilian Actors
Plot conflicts involving civilians
codes <- c('17', '27', '37', '47', '78')
par(mfrow = c(3,2), mar = c(0,0,1,0)) 

for(code in codes) {
  # Filter for the specific interaction type
  filtered_data <- conflict_data_sf %>%
    filter(interaction_code == code)
  
  # Convert to ppp object with boundary window
  boundary_owin <- as.owin(boundary_sf)
  ppp_obj <- as.ppp(filtered_data$geometry, W = boundary_owin)
  
  # Rescale the ppp object
  ppp_obj_km <- rescale(ppp_obj, 1000, "km")
  
  # Compute KDE
  kde <- density(ppp_obj_km,
                 sigma = 71.831,
                 edge = TRUE,
                 kernel = "quartic")
  
  # Plot the KDE
  plot(kde, main = paste(interaction_map[code]), col = colours)
}

Observations

Across all actors, we see that perhaps civilians are involved with the biggest range of actors.

Notably, we see the lowest intensity of clashes with identity militia which could stem from how civilians might be part of the same social/ethnic networks as the identity militia and these militia may seek the support of civilian populations too.

In contrast, state forces and civilians have the highest intensity of armed conflicts

  • In Western Rakhine State, Buddhist-majority and Muslim-minority communities have experienced significant tension causing state forces to act in response to these tensions which has led to severe clashes with civilian populations.

  • Reports of human rights abuses, including arbitrary arrests, torture, and extrajudicial killings by state forces, has led to large-scale displacement and humanitarian crises, further exacerbating tensions between state forces and civilians. The struggle for resources and safety often intensifies the conflicts. (Farge E. & Mantovani C., 2024)

4. Civilians Involvement by Event Type
Compute KDE of Civilians by Event Type
# Set Up
par(mfrow = c(2,2), mar = c(0,0,1,0)) 
codes <- c('17', '27', '37', '47', '78')

conflict_data_sf %>%
  filter(interaction_code %in% codes) %>%
  group_by(event_type) %>%
  group_split() -> conflict_by_event_type

# Convert the sf object to owin
district_boundary <- as.owin(st_as_sfc(boundary_sf))

kde_list <- lapply(seq_along(conflict_by_event_type), function(i) {
  data <- conflict_by_event_type[[i]]
  event_type <- unique(data$event_type)
  ppp_obj <- as.ppp(st_geometry(data), W = district_boundary)
  ppp_obj <- rescale(ppp_obj, 1000, "km")
  kde <- density(ppp_obj,
                 sigma=71.831,
                 edge=TRUE,
                 kernel="quartic")
  plot(kde, main = paste(event_type), col=colours)
  return(kde)
})

Observations

Civilians are seen to be embroiled mostly in conflicts resulting from strategic development and violence against civilians events, primarily in Central and Western Myanmar. Fortunately, conflicts involving explosions or remote violence is not as intense against civilians but this raises a great concern on the humanitarian crisis faced by civilians in Myanmar.

5.5.5 OpenStreetMap Myanmar - Spatial Points

Using tmap functions, I will display an interactive view of te KDE layers on openstreetmap of Myanmar to observe intensity of conflicts at the district level. We can see higher intensities of conflict in central districts (e.g. Yinmarbin, Pakokku and Shwebo) and western districts (e.g. Yangon).

ppp_obj<- as.ppp(conflict_data_sf$geometry)
ppp_owin_obj <- ppp_obj[myanmar_owin]
ppp_owin_obj_km <- rescale(ppp_owin_obj, 1000, "km")

kde_fixed <- density(ppp_owin_obj_km, sigma=bw.CvL, edge=TRUE, kernel="quartic")

raster_kde_fixed <- raster(kde_fixed)

projection(raster_kde_fixed) <- CRS("+init=EPSG:32647 +units=km")
# Plot KDE Map on OpenStreetMap
tmap_mode('view')
tmap mode set to interactive viewing
kde_fixed_output <- tm_basemap(server = "OpenStreetMap.HOT") +
  tm_basemap(server = "Esri.WorldImagery") +
  tm_shape(raster_kde_fixed) +
  tm_raster("layer",
            n = 10,
            title = "KDE Fixed CvL",
            alpha = 0.6,
            palette = c("#fafac3","#fd953b","#f02a75","#b62385","#021c9e")) +
  tm_shape(boundary_sf)+
  tm_polygons(alpha=0.1,id="DT")+
  tmap_options(check.and.fix = TRUE)

kde_fixed_output
tmap_mode("plot")
tmap mode set to plotting

6. 2nd Order Spatial Point Patterns Analysis

Unlike 1st-order analysis, which studies the intensity of points (e.g., density), let’s also leverage 2nd-order analysis to examine how points are distributed relative to each other, which can offer deeper insights into the spatial interaction between events.

I will use the K-function and L-function is to understand the spatial relationships between events, particularly focusing on whether the points exhibit clustering, uniformity, or randomness.

6.1 Using K-Function Estimation

K-function helps detect spatial patterns by comparing the observed distribution of points against a random pattern at different distances.

6.1.1 Yinmarbin District

1) Computing K-Function Estimation

For Yinmarbin district, let’s compute K-function estimates by using Kest() of the spatstat package by using the yinmarbin_ppp_owin masked ppp object we created previously.

K_ck = Kest(yinmarbin_ppp_owin, correction = "Ripley")
plot(K_ck, . -r ~ r, ylab= "K(d)-r", xlab = "d(km)",
       main = paste("Yinmarbin District (K-Function)"))

Observations

How to interpret the plot:

  • K-iso represents the observed or estimated K-function value calculated from the actual data
  • K-pois is the theoretical K-function that represents the expected K-function

With that said…

We can observe how the observed line (K-iso) is found to lie above the theoretical line (K-pois) which suggests conflict points in Yinmarbin are highly clustered. Clustering is the strongest at a 20km distance which suggests large-scale clustering. In fact, it is more clustered together than expected by the null hypothesis.

Note: Since I had used the default edge = TRUE settings, edge correction will account for missing neighbours outside the boundary which helps maintain an accurate estimate of the K-function. Hence, there is marginal difference in the actual and expected K-function.

2) Performing Complete Spatial Randomness Test

To confirm the observed spatial patterns above, a hypothesis test (i.e. Monte Carlo simulation test) will be conducted. The hypothesis and test are as follows:

  • Ho = The distribution of conflict events in Myanmar are randomly distributed.
  • H1= The distribution of conflict events in Myanmar are not randomly distributed.
  • The null hypothesis will be rejected if the observed K-function lies above/below the theoretical K-function and envelope.

By using envelope(), we can get a more robust interpretation by comparing the observed K-function against a simulation envelope of K-functions generated under the null hypothesis.

Note

To achieve a 95% confidence envelope in a K-function test with Complete Spatial Randomness, I will need to exclude the upper 2.5% and lower 2.5% of the simulated K-functions., i.e. I will need to generate at least 40 simulations where nsim = 39.

We’ll set rank = 1 for a conservative setting, such that the envelope is based on the extreme values (highest and lowest) from the simulations.

# Monte Carlo test with K-function
K_yinmarbin_csr <- envelope(yinmarbin_ppp_owin, Kest, 
                            nsim = 39, rank = 1, glocal=TRUE)
plot(K_yinmarbin_csr, main = paste("Yinmarbin District (CSR)"))

Observations

The tight envelope suggests we can be confident that deviations between the observed and theoretical lines are meaningful, rather than due to random variation. It also indicates that simulations under CSR are producing homogeneous patterns.

On the other hand, the observed K-function line is constantly above the shaded region of the envelope and theoretical line, suggesting that observed spatial pattern are more clustered than expected under the null hypothesis of CSR

Just to be sure that there is no improvement in a higher number of simulations, I ran it again with nsim = 99. As shown below, the k-function output is similar but with a somewhat smoother plot than before since we are working with less variability. Hence, I’ll stick to nsim = 39.

6.1.2 Shwebo District

1) Computing K-Function Estimation

For Shwebo district, let’s compute K-function estimates by using Kest() of the spatstat package.

K_ck = Kest(shwebo_ppp_owin, correction = "Ripley")
plot(K_ck, . -r ~ r, ylab= "K(d)-r", xlab = "d(km)",
       main = paste("Shwebo District (K-Function)"))

The observed line is consistently above at all distances, implying that conflict events are clustered at both small and large areas, rather than being randomly distributed across the area.

Observations

The observed line is consistently above at all distances, implying that conflict events are clustered at both small and large areas, rather than being randomly distributed across the area.

2) Performing Complete Spatial Randomness Test

# Monte Carlo test with K-function
K_shwebo_csr <- envelope(shwebo_ppp_owin, Kest, 
                     nsim = 39, rank = 1, glocal=TRUE)
plot(K_shwebo_csr, main = paste("Shwebo District (CSR)"))

Observations

Since the observed line lies above the envelope and the shaded region of the envelope is rather narrow, it strongly suggests clustering in Shwebo, rather than due to random variation.

6.1.3 Pakokku District

1) Computing K-Function Estimation

For Pakokku district, let’s compute K-function estimates by using Kest() of the spatstat package.

K_ck = Kest(pakokku_ppp_owin, correction = "Ripley")
plot(K_ck, . -r ~ r, ylab= "K(d)-r", xlab = "d(km)",
       main = paste("Pakokku District (K-Function)"))

Observations

We can observe how the observed line (K-iso) is found to be above the theoretical line (K-pois) which suggests conflict points in Pakokku are highly clustered. In fact, it is more clustered together than expected by the null hypothesis.

2) Performing Complete Spatial Randomness Test

# Monte Carlo test with K-function
K_pakokku_csr <- envelope(pakokku_ppp_owin, Kest, 
                     nsim = 39, rank = 1, glocal=TRUE)
plot(K_pakokku_csr, main = paste("Pakokku District (CSR)"))

Observations

Likewise, we can be certain about the clustering patterns in Pakokku since the observed K values deviates above the envelope as shown.

6.1.4 Mandalay District

1) Computing K-Function Estimation

For Mandalay district, let’s compute K-function estimates by using Kest() of the spatstat package.

K_ck = Kest(mandalay_ppp_owin, correction = "Ripley")
plot(K_ck, . -r ~ r, ylab= "K(d)-r", xlab = "d(km)",
       main = paste("Mandalay District (K-Function)"))

Observations

Significant clustering is observed generally across multiple distances in Mandalay and clusterinng is particularly evident at longer spatial distances in Mandalay seen from how the deviations of observed K values increase by the distance.

2) Performing Complete Spatial Randomness Test

# Monte Carlo test with K-function
K_mandalay_csr <- envelope(mandalay_ppp_owin, Kest, 
                     nsim = 39, rank = 1, glocal=TRUE)
plot(K_mandalay_csr, main = paste("Mandalay District (CSR)"))

Observations

The envelope size continues to be narrow with K values plotted above the envelope. This strongly suggests that the point pattern is not random, and the points are clustered.

6.2 Using L-Function Estimation

In this section, I will be computing L-function via Lest() of spatstat package which is normalises the K-function to a linear scale for easier interpretation.

6.2.1 Yinmarbin District

1) Computing L-function Estimation

L_yinmarbin = Lest(yinmarbin_ppp_owin, correction = "Ripley")
plot(L_yinmarbin, . -r ~ r, ylab= "L(d)-r", xlab = "d(km)",
     main = paste("Yinmarbin District (L-Function)"))

Observations

The theoretical L-function under CSR remains at 0, meaning that if the points were randomly distributed, the L-function would be expected to show no significant clustering or dispersion.

There is a consistent deviation of each observed L value from the theoretical value across all distances, particularly with the strongest deviation at <5 km and > 20km. This confirms that conflict points are both localised in towns of Yinmarbin and found across large areas in Yinmarbin.

2) Performing Complete Spatial Randomness Test

To confirm the observed spatial patterns above, a hypothesis test will be conducted. The hypothesis and test are as follows:

  • Ho = The distribution of conflict events in Myanmar are randomly distributed.
  • H1= The distribution of conflict events in Myanmar are not randomly distributed.
  • The null hypothesis will be rejected if the observed L-function lies above/below the theoretical L-function and envelope.

I will also perform monta carlo simulation test using envelope() of the spatstat package.

# Monte Carlo test with L-function
L_yinmarbin_csr <- envelope(yinmarbin_ppp_owin, Lest, 
                     nsim = 39, rank = 1, glocal=TRUE)
plot(L_yinmarbin_csr, . - r ~ r, xlab="d(km)", ylab="L(d)-r",
       main = paste("Yinmarbin District (CSR)"))

Observations

The envelope outputted is narrow which once again reflects lower variability in spatial distribution in Yinmarbin where a large number of points found in this district could have led to more stable results. There is also statistically significant clustering particularly in specific towns in Yinmarbin and across Yinmarbin itself, since we see larger deviations of observed L values from the envelope.

6.2.2 Shwebo District

1) Computing L-function Estimation

L_shwebo = Lest(shwebo_ppp_owin, correction = "Ripley")
plot(L_shwebo, . -r ~ r, ylab= "L(d)-r", xlab = "d(km)",
     main = paste("Shwebo District (L-Function)"))

Observations

The observed L values obtained for Shwebo districct tend to show a smaller deviation above the theoretical L-function in smaller areas like towns, suggesting localised clustering is less pronounced. Conversely, larger deviations are seen at larger distances which shows broader patterns of clustering are present in Shwebo district.

2) Performing Complete Spatial Randomness Test

# Monte Carlo test with L-function
L_shwebo_csr <- envelope(shwebo_ppp_owin, Lest, 
                     nsim = 39, rank = 1, glocal=TRUE)
plot(L_shwebo_csr, . - r ~ r, xlab="d(km)", ylab="L(d)-r",
       main = paste("Shwebo District (CSR)"))

Observations

The envelope size continues to grow larger as distance increases, which indicates greater variability in spatial patterns of the Shwebo district but overall, we can confirm strong clustering of conflicts especially in bigger areas of the district which suggests wide-scale hotspots, as well as some localised clustering at smaller scales.

6.2.3 Pakokku District

1) Computing L-function Estimation

L_pakokku = Lest(pakokku_ppp_owin, correction = "Ripley")
plot(L_yinmarbin, . -r ~ r, ylab= "L(d)-r", xlab = "d(km)",
     main = paste("Pakokku District (L-Function)"))

Observations

Like the other districts, there is clustering seen in Pakokku across smaller areas like towns and the distribution across the district itself.

2) Performing Complete Spatial Randomness Test

# Monte Carlo test with L-function
L_pakokku_csr <- envelope(pakokku_ppp_owin, Lest, 
                     nsim = 39, rank = 1, glocal=TRUE)

Observations

Under the CSR test, the intensity of clustering weakens across larger distances (>20km) than we would have assumed in the previous step. This happens as Monte Carlo simulation accounts for random fluctuations especially at larger distances (e.g., >20 km) where random variability in the point pattern naturally increases.

In short, significant clustering is present at both small and large scales, but with diminishing intensity as we measure larger distances.

6.2.4 Mandalay District

1) Computing L-function Estimation

L_mandalay = Lest(mandalay_ppp_owin, correction = "Ripley")
plot(L_mandalay, . -r ~ r, ylab= "L(d)-r", xlab = "d(km)",
     main = paste("Mandalay District (L-Function)"))

Observations

The distance scale for Mandalay is the smallest (0-10km) since the district itself has a smaller area. Nonetheless, we see almost a consistent spread of localised and widespread distribution of conflict points in Mandalay which indicates that conflicts occur across both specific towns and the across the districct.

2) Performing Complete Spatial Randomness Test

# Monte Carlo test with L-function
L_mandalay_csr <- envelope(mandalay_ppp_owin, Lest, 
                     nsim = 39, rank = 1, glocal=TRUE)
plot(L_mandalay_csr, . - r ~ r, xlab="d(km)", ylab="L(d)-r",
       main = paste("Mandalay District (CSR)"))

Observations

The Monte Carlo simulation shows that there’s statistically significant clustering across all distances since the observed L values consistently lie above the envelope. There is also a marginal increase in variability of random points observed which is inherent in larger distance scales as seen from the slight widening of the envelope as distance increases.

Continue to Part 2 >